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Abstract 



We study the relationship between long cycles and Bose-Einstein condensation (BEC) in 
the case of several models. A convenient expression for the density of particles on cycles of 
length q is obtained, in terms of q unsymmetrised particles coupled with a boson field. 

Using this formulation we reproduce known results on the Ideal Bose Gas, Mean-Field and 
Perturbed Mean-Field Models, where the condensate density exactly equals the long cycle 
density. 

Then we consider the Infinite-Range-Hopping Bose-Hubbard Model: 
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in two cases, first for A = +00, otherwise known as the hard-core boson model; and secondly 
for A finite, representing a finite on-site repulsion interaction. 

For the hard-core case, we find we may disregard the hopping contribution of the q unsym- 
metrised particles, allowing us to calculate an exact expression for the density of particles 
on long cycles. It is shown that only the cycle of length one contributes to the cycle density. 
We conclude that while the existence of a non-zero long cycle density coincides with the 
occurrence of Bose-Einstein condensation, the respective densities are not necessarily equal. 

For the case of a finite on-site repulsion, we obtain an expression for the cycle density 
involving the partition function for a Bose-Hubbard Hamiltonian with a single-site correction 
again by neglecting the q unsymmetrised hop. Inspired by the Approximating Hamiltonian 
method we conjecture a simplified expression for the short cycle density as a ratio of single- 
site partition functions. In the absence of condensation we prove that this simplification is 
exact and use it to show that in this case the long-cycle density vanishes. In the presence 
of condensation we can justify this simplification when a gauge-symmetry breaking term is 
introduced in the Hamiltonian. Assuming our conjecture is correct, we compare numerically 
the long-cycle density with the condensate and again find that though they coexist, in general 
they are not equal. 
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Chapter 1 
Introduction 



Summary 

To begin it is useful to look back through the literature to trace the progression of discoveries 
made by many great minds on the phenomenon of Bose-Einstein condensation in quantum 
statistical mechanics. We discuss the order parameters that have been postulated to indicate 
the presence of this condensate, and give details on all known about their relationships to 
date. Finally we briefly detail the layout of this thesis. 



1.1 Bose-Einstein Condensation 

In 1925, the phenomenon of Bose-Einstein condensation was predicted by Einstein following 
correspondence with S.N. Bose, who had found a novel derivation of photon statistics and 
the Planck distribution [36]. Einstein extended these ideas to an ideal monatomic gas [37] 
and observed that below some critical temperature, a number of molecules which grows with 
the total density, transition to an "unmoving state," and this condensed substance coexists 
with the monatomic gas. These papers appeared about a year before the development of 
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Quantum Mechanics and was the first time the wave nature of matter was explored, with 
particles obeying the same statistics as photons. 

Since then we now know that all particles in nature can be classified into two classes according 
to their intrinsic angular momentum or spin. Those with half-integer spin are called fermions, 
those with integral spin are bosons. Fermions are restricted by the Pauli- Exclusion principle, 
which prevent multiple fermions from existing in the same quantum state, however bosons 
have no such restriction. The aforementioned "unmoving state" is now known to consist 
of bosons at (or near) the bosonic ground state, with particles in this state denoted as the 
(generalised) Bose-Einstein condensate (BEC). 

Nevertheless at the time, Einstein's claim of a phase transition in the ideal bose gas was 
dismissed, largely due to the criticism in 1927 of G. Uhlenbeck [83] (a graduate student of 
P. Ehrenfest) who argued that a gas of finitely many particles can exhibit no discontinuity, 
which was technically correct at the time, as it was yet to be shown that singularities in 
thermodynamic functions are a consequence of the thermodynamic limit (i.e. second-order 
phase transitions). This opinion prevailed until 1937 when H.A. Kramers pointed out the 
importance of the limit, causing Uhlenbeck [84] to withdraw his objections and assert that 
Einstein's formula for the particle number-density is correct in the thermodynamic limit. 

Meanwhile although the experimentalist H. Kamerlingh Onnes [85] had succeeded in liqui- 
fying Helium ( 4 He) in 1908, it was not until 1928 that his student W.H. Keesom with M. 
Wolfke [55] discovered that at temperatures below 2.18K, a second phase of liquid helium 
appears, and christened it Helium II. Four years later Keesom and K. Clusius [54] measured 
the famous peak in the specific heat of helium, whose shape inspired the name A-transition 
for the behaviour at the critical temperature T\ = 2.18K (see Figure 1.1). Further research 
[53; 3] in 1938 into Helium II found that a portion of the liquid flows without any apparent 
viscosity, a phenomenon now called superfluidity. In 1937 F. London [62] connected this 
curious phase transition in 4 He (a bosonic gas) with the disputed ideas of Einstein on the 
ideal gas. He noted that Einstein's formula for the transition temperature was a good estim- 
ate for T\ and consequently introduced the concept of macroscopic occupation of the ground 
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Figure 1.1: The X-transition in 4 He. 

state and related it to the coherence properties of the condensate. That the superfluidity 
of 4 He is due to the presence of Bose-Einstein condensation is now a fact which is largely 
accepted today. 

Cowley and Woods [23] in 1968 were the first to experimentally detect condensate at T = 
1.1K. Later work [2; 27] repeated and improved these experiments and found the critical 
temperature to be 2.24 ± 0.04K, closely matching the corresponding value for the lambda- 
point T\, lending weight to a correspondence between condensation and superfluidity. A 
more thorough historical account of these discoveries can be found in [47]. While BEC 
has been proved experimentally, nevertheless it remains a valuable concept to the quantum 
statistical mechanics community as it is the only known phase-transition taking place in the 
absence of interactions, its origin being purely quantum mechanical. 

A general theory for condensation of interacting bosonic gasses remain to this day a no- 
toriously difficult problem. The theory of the ideal gas has led to the formulation that 
"Bose-Einstein condensate" labels those particles which occupy lowest kinetic energy state, 
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A rigourous proof that Qq > for systems with realistic inter-particle interactions is an 
ongoing pursuit. In 1947 N.N. Bogoliubov computed the excitation spectrum for the weakly 
interacting Bose gas and pointed out the similarities with the spectrum of Helium II. This 
implied that the inclusion of interactions in the model had little effect upon the condensate, 
but instead had a pronounced impact upon the wave function. 

Until recently in fact, BEC had not been proved for many-body Hamiltonians with genuine 
interactions except for one special case: hard core bosons on a lattice at half-filling, i.e., 
number of particles is exactly half the number of sites, in 1978 [35; 56]. Five years later in 
1983, E. Buffet, Ph. de Smedt and J.V. Pule [20] showed the existence of BEC for a class 
of homogeneous bosonic systems with interactions close to the van der Walls limit with a 
gap in the one-particle excitation spectrum. Independently in 2002, E. Lieb and R. Seiringer 
[60] considered a model of an interacting system of a trapped (and thus inhomogeneous) 
gas in the Gross-Pitaevski limit. J. Lauwers, A. Verbeure and V.A. Zagrebnov [58] in 2003, 
proved the occurrence of zero-mode Bose-Einstein condensation for a class of continuous 
homogeneous systems of boson particles with super-stable interactions and a one-particle 
energy gap. 

Because of the difficulty of treating the boson gas interacting through a pair potential, 
various approximations have been proposed. One approach is to isolate the portion of the 
interaction which depends only upon the occupation numbers, separating the interaction 
it into a diagonal part (with respect to the occupation number basis) and an off-diagonal 
part (see [59]). K. Huang, C.N. Yang and J.M. Luttinger in 1957 [50] argued that some of 
the diagonal terms can enhance the amount of condensate, while Bogoliubov and his school 
claimed that it contains terms which deplete it. C.N Yang and CP. Yang [98] considered 
a special case in d = 1 of a repulsive delta potential and prove that it exhibits no phase 
transition (see also [32; 31]). However considering the diagonal portion of this interaction 
alone, the Hamiltonian reduces to that of the Huang- Yang-Luttinger model which has a 
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guaranteed phase transition at d — 1. This points to the idea that the off-diagonal terms 
must suppress this condensate, which indicates a more complex interplay of the diagonal and 
off-diagonal portions of the interaction than previously thought. Select models considering 
only the diagonal terms include the Mean-Field model [24; 91], the Huang- Yang-Luttinger 
model [89; 88], the Perturbed Mean-Field model [87; 28], and the Full-Diagonal model [29], 
all of which can possess a critical phase. 

Another model which has been considered includes non-diagonal BCS-type interaction terms 
in addition to the diagonal terms [70; 100]. 

Analysis of further models was made possible by two Ansatze by Bogoliubov, the first to 
truncate the full Hamiltonian of the interacting bosons to produce the Weakly Imperfect 
Bose Gas (WIBG)[99], the second to perform apt substitutions of some operators by c- 
numbers, known as the Bogoliubov approximation. Condensation is shown to appear in 
the WIBG (with some peculiarities). Investigations into the experimentally observed Mott 
insulator- superfluid (or condensate) phase transition [44] resulted in the analysis of the so- 
called Bose-Hubbard model, proposed by Fisher et. al. [43], and successfully analysed using 
the Bogoliubov approximation in [17; 30] to exhibit BEG 

1.2 The Feynman-Kac Integral and Cycle Lengths 

An alternative parameterisation to predict the occurrence of the condensate was postulated 
by R.P. Feynman in 1953 [42]. He was discouraged by the fact that energies of a complex sys- 
tem are "hopelessly" difficult to calculate, so he tried combine his path integral formulation 
with Bose statistics to see if a transition could be found by other means. 

The idea behind the Feynman path integral can be traced back to a paper by Dirac [26] in 
1933 on an alternative formulation of quantum mechanics. Initially the quantum theories 
concurrently developed by Heisenberg and Schrodinger and others in the 1920 were built as 
a close analogy to classical mechanics where they considered the Hamiltonian H, an energy 
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function of coordinates and momenta. Dirac found that the Lagrangian L, a function of the 
coordinates x(t) and velocities x(t), to be more intuitive than the Hamiltonian approach. 
His reasoning is as follows, in classical dynamics the action is given by 

S[x(-)} = J L(x(t),x(t))dt (1.1) 

where the integral is taken between initial and final times and t. The action must depend on 
the path x( ■ ) taken by the particles, and hence is a functional of these paths. The principle 
of least action states that for small variations of these paths, the action S is minimised. 
He hoped that these classical paths could somehow be useful to describe the evolution of a 
quantum particle. 

In his article, Dirac considers the time propagation operator Q- ltH / h , For two points x and y 
in a region A, quantum mechanics tells us that we can write the probability P t (x,y) that a 
particle with initial position x propagates to position y after time t, as Pt(x, y) = \K t (x, y)\ 2 , 
the square of the amplitudes 

K(x,y;t) = (y\e- UH / h \x). 

This amplitude is the sum of the amplitudes of all possible paths which start at x and end 
at y after time t. We will denote such a path by cu* : [0,t] — > A and denote its amplitude 
by lf[c<;* ]. Let fi* be the set of all these paths. Then 

K(x,y;t)= £ K[ulJ. 

We would like to see how each path contributes to the total amplitude. Dirac suggests 
that the amplitude -^[cu* ] "corresponds to" the quantity 

exp {^[<^)]} 
and later writes that it is the "classical analogue" of the propagator. 
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By 1948 Feynman [41] had developed this suggestion and showed it was more than an 
analogy: it was an equality. He succeeded in formulating his "space-time approach," deriving 
a third formulation of Quantum Mechanics using this expression for the propagator in terms 
of the action (a formulation later shown to be equivalent to the Heisenberg and Schrodinger 
pictures). And while Dirac considered only the classical path, Feynman showed that all 
paths contribute: in a sense, the quantum particle takes all paths, and the amplitudes for 
each path add according to the usual quantum mechanical rule for combining amplitudes. 
What was lacking was a rigourous definition of summing over these paths. 

Surprisingly a well-defined calculus of such paths was already well known to mathematicians 
by the 1920s, due to N. Wiener and his study of stochastic processes and Brownian motion 
(extending the earlier work of Einstein and M. Smoluchowski). He introduced the concept of 
integration in functional spaces that now bears his name: the Wiener measure and Wiener 
integral. In 1949, Kac [51] considered Feynman's path integral and realised that if one 
replaces the time parameter in the Dirac- Feynman expression for the action, equation (1.1), 
with a purely imaginary time, the Feynman path integral becomes a Wiener integral which 
(with some work) can be written in terms of the Hamiltonian H = H + V: 



K(x,y-t^-ihf3) = (y\e^ H \x) = J Wg y (du) exp j- jf 





V(u(s))ds 



which is now called the Feynman-Kac path-integral. 

Returning to Feynman's proposition of an order parameter for BEC, he started out by 
analysing the partition function for a system of n particles interacting via a pair potential 
of the form V(x,y) = V(\x — y\). Using the Feynman-Kac formula, one can write the e~ /3H , 



K(x u ... J Wl^idLJi) ■ • ■ J W x l yn (du n ) 

x exp j- Y^J* Vfais) - w fc (s))ds j (1.2) 

where the Wiener measure W£ y (uj) describes a Brownian path that starts at x and ends at 
y with "time" (3, i.e. lo(0) = x and uj(/3) = y. 
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The Feynman path integral has become a very powerful tool in Statistical Mechanics. For a 
system of bosons with Hamiltonian H, the canonical partition function can be written as 

Z v (n) = trace [o+e _/3H ] = ^ J dxi ' ■■ dx nK(x 1 , . . . , x n ; x^, . . .,x ff ( n )) (1.3) 

where the sum over permutations appears because we take the trace over the symmetric 
states (as indicated by cr + ), and we integrate over all particle positions. (See Appendix A of 
[81] for a clear derivation of this expression). 

P 




y 



Figure 1.2: An interpretation of the Feynman-Kac representation of the partition function 
for a gas of bosons. The horizontal plane represents the d spatial dimensions, and the vertical 
axis is the imaginary time direction. This picture considers five particles on two cycles, one 
of length 1 and the other of length 4 (Picture from [81]). 

In Statistical Mechanics, /3 is the inverse temperature. However if one were consider it as a 
time parameter for the Feynman-Kac formula, then one returns to the intuitive picture of 
moving particles following paths. Then equation (1.2) can be thought of as n particles with 
initial positions Xi and final positions ?/, at "time" (3. 

Now if one considers the integrand of the partition function expression (1.3), one sees that 
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the initial and final positions for the particles are fixed, but a particle with initial position 
Xi will have final position x n ^, thus its path will depend upon the permutation n. 

Figure 1.2 provides an example, where five particles propagate in the imaginary time direction 
on paths with endpoints specified by the permutation 



The particle labelled five is on a closed Brownian bridge, starting and ending at the point 
£5. We shall say this is a cycle of length one, whereby it took one timestep /3 for the particle 
to return to its original position. 

However the particle labelled one however lands on a different site, X2, after time /3. Should 
it continue its journey, it will land at £4 at time 2/3, x 3 at time 3/3 and back to its starting 
position x\ after four timesteps. We shall call this a cycle of length four. 

Feynman considered the probability that a given particle belongs to a cycle of length, say 
q, and questioned if there was a strictly positive probability that, after taking the thermo- 
dynamic limit, a cycle of infinite length exist? He postulated that the occurrence of infinite 
(long) cycles is indicative of the presence of BEC. 

Three years after Feynman's conjecture, O. Penrose and L. Onsager [67] introduced a different 
order parameter for BEC, related to the concept of "off-diagonal long-range order" (OLDRO) 
[97] which investigates the correlation between two particles at positions x and y, i.e. 



We say the system displays OLDRO when D p (x, y) > when \x — y\ —> 00. This correlation 
can be expressed in the Feynman-Kac representation, which involves an open cycle starting 
at x and finishing at y which may wind many times around the imaginary time direction 




4 1 3 



2 3 4 
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(see Figure 1.3): 

D p (x,y) = lim - } — -. — 1 —— [ dx 2 . . . dx n V / dW? * M ■ ■ ■ / dW? * (w n ) 
pv yy »,v-*» Z v (n)V (n - 1)! J ^ J x ^m y u J v n ^ 

Tl/V=p TTgSn 

x exp J - V /" ^(^(s) -w fc (s))dsl (1.4) 
I J 

where we have set X\ — x,Xi — y and Xj = Xj for j — 2, . . . ,n. As |x — y\ — > oo there 
corresponds a notion of infinite winding which is reminiscent of Feynman's approach. Penrose 
and Onsager observed that there should be BEC when the fraction of the total particle 
number belonging to infinite cycles is strictly positive. 

(3 



y 




y 



Figure 1.3: Interpretation of Off-Diagonal Long-Range Order in terms of Feynman-Kac 
integrals 

It was generally agreed at the time that this parameter is the correct one for BEC, however 
infinite cycles were still considered noteworthy, with for instance a possible role in super- 
fluidity [68], and showing analogous phase transitions to various physical processes such as 
percolation, gelation and polymerisation, see for example [22; 75; 76; 25]. 
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However it was not until 1993 that A. Siito [78; 79] gave a precise mathematical and quantit- 
ative formulation relating the occurrence of BEC (and hence OLDRO) with the existence of 
infinite cycles. He proved in the cases of the free bose gas and the mean field gas that in fact 
infinite cycles do occur below the transition temperature for BEC. In 2001, D. Ueltschi [82] 
considered a simple lattice model in which the integral over Brownian paths is replaced by 
an effective weight. He proved that infinite cycles do not occur for high temperatures. Later 
in 2003, T.C. Dorlas, Ph. Martin and J.V. Pule [34] considered a model of a Bose gas with 
some genuine short range interactions between modes using the theory of Large Deviations 
developed by Varadhan [95], and proved that the average fraction of infinite cycles is equal 
to the density of condensate. 

However more recent investigations cast doubt upon the suitability of infinite cycles as an 
order parameter for BEC. Ueltschi [82] argues that both order parameters are identical in the 
case of a weakly interacting Bose gas, but in the presence of a strong inter-particle interaction 
this correlation breaks down. The example offered is a crystalline system. The off-diagonal 
correlation function decays exponentially as the inter-particle distance increases so there 
in no OLDRO. However Brownian motion occurs as a consequence of the interchange of 
particles via quantum tunnelling. Due to the strong interaction, these Brownian paths need 
to avoid each other, which for dimensions greater than two and for low enough temperatures 
can result in a macroscopic number of particles on random walks of infinite length. 

As a result there exists uncertainty as to exactly what connection there is between Bose- 
Einstein condensation and infinite cycles, which motivates the subject of this thesis. Here 
we compare the long cycle density with that of the condensate for two models of interacting 
bosons, both of which involve an on-site repulsive interaction. 

We find a weaker result than [34], that the existence of long cycles coincides with the presence 
of BEC, but their respective densities are not necessarily equal. In fact, the long cycle density 
can be both greater than or less than the condensate density. It appears that the single-site 
repulsion term can either increase or decrease the quantity of long cycles, depending upon 
the thermodynamic parameters of the system. Nevertheless this analysis does lend weight 
to Feynman's conjecture in the absence of strong interactions. 
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Outline of Thesis 

The object of this thesis is to study the relation between the occurrence of BEC and the 
existence of infinite cycles for some Boson lattice models, in particular the infinite-range- 
hopping Bose-Hubbbard model with both hard cores and without. These models correspond 
to bosonic systems that describe particles hopping on a lattice and interacting among them- 
selves with the constraint that their wavefunctions are invariant under permutations of the 
particles. The partition function can induce a probability measure upon the set of permuta- 
tions, with which one can calculate expectations of cycle lengths. 

In Chapter 2 we apply the techniques for cycle statistics described in [34], following [64], to 
establish the general framework for evaluating the density of particles on cycles of particular 
lengths, and show that this formulation uniquely partitions the particles into their cycle 
structures. 

In Chapter 3 we formulate a relation between occupation numbers and the cycle structure 
of models which satisfy particular conditions, which implies the condensate density equals 
the infinite cycle density. We then check the conditions in the case of three simple models, 
the Ideal Bose Gas, the Mean-Field Bose Gas and the Perturbed Mean-Field Bose gas, to 
reproduce the results of Siito [78; 79] and Dorlas, Martin and Pule [34]. 

In Chapter 4 we consider the Infinite-Range-Hopping Bose-Hubbard Model with hard-cores, 
i.e. infinite on-site repulsion to prevent multiple particles occupying individual sites. We 
calculate the density of particles on cycles on infinite length and compare the result with the 
condensate density. 

In Chapter 5 we consider the Infinite-Range-Hopping Bose-Hubbard Model (without hard- 
cores) where now there is a finite on-site repulsion to discourage multiple particles per site. 
Again we consider infinite cycle density and compare the result with the condensate density. 
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Summary 

In this chapter we will introduce a general framework for cycle statistics by expanding the 
bosonic partition function into a sum over permutations of an unsymmetrised trace term and 
inducing a probability measure upon these permutations. By decomposing these permutations 
into their cycle structures, we derive a convenient expression for the density of particles on 
cycles of particular lengths, and show that this formulation uniquely partitions the particles 
depending upon the cycle lengths they live on. In taking the thermodynamic limit we define 
the long cycles to be those which are of infinite length, and short cycles as those which are 
finite in length. We find that we can uniquely decompose the total density of particles into 
those on long cycles and those on short cycles, i.e. p = fthort + £?iong- Finally we also obtain 
an expression for the expectation of certain operators in terms of cycle statistics, which will 
be useful for analysing occupation number densities and off-diagonal long-range order in later 
chapters. 
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2.1 Cycle Representation of the Partition Function 

In his paper investigating the low-temperature properties of liquid He 4 , R. Feynman [42] 
re-expressed the partition function of a boson gas in terms of the cycle structures of per- 
mutations induced by Bose-Statistics. He postulated that at low temperatures (or high 
densities), bosons belong to cycles of different lengths, and below some critical temperature 
cycles of infinite length appear. Feynman associated this critical behaviour with the phase 
transition of He 4 known as Bose-Einstein condensation. 

In this chapter we will establish a general framework for cycle statistics as described in [34], 
following [64]. This framework will be applied throughout this script so it is valuable to 
lay down the notational foundation for this thesis before applying the techniques to various 
models. 

To begin, we shall consider a system of n identical bosons enclosed in a region A with volume 
|A|, in thermal equilibrium at temperature T (i.e. in the canonical ensemble). Denote the 
inverse temperature (5 := ksT, where ks is Boltzmann's constant. For this ensemble we take 
the thermodynamic limit at a fixed temperature f3 by letting the quantities n and |A| tend 
to infinity in such a way that n/] A] = p for a fixed p. For the entirety of this text we shall 
fix (3, and drop it from our notation where possible. 

Denote the single particle Hilbert space by T-L h . Let / be the identity operator on T-L A . For an 
operator A on H A denote as the operator on an unsymmetrised (i.e. no Bose-statistics) 
n-particle Hilbert space 



When considering bosons, we need to restrict the system to the symmetric subspace of 7i A . 
Let S n be the set of all permutations of {!,... , n}. The symmetrisation projection a± on 



n times 



A {n) = A®I®---(g)I + I(g)A(g)---®I-\ h / <8) • • • <8> / <8> A. 



% A is defined by 
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where U n : U { ? ] H- is a unitary representation of the permutation group S n on "H 
defined by 

[4(01 ® 02 <E> • • • ® 4>n) = 0tt(1) ® 4> n (2) ® • • • <S> 0,r(n) 



(n) 
A 



with 



G 7/ A , J = 1, • • - ,n and any tt G S^. Thus we denote the symmetric or bosonic 



n-particle Hilbert space by % 



in) 



°H A n) . We assume that there is a symmetric Hamilto- 



nian H2 on T-L^ whose restriction to "% A n + is the Boson Hamiltonian for our system. The 
canonical partition function for the n-particle Hamiltonian if™ can be written as 



ZAn) = trace 

rt A,+ 



trace 

rt A 



ly 

7T6S n 



trace 

ft A 



This leads us naturally to define a probability measure on S n by 



£4e' 



-PHI 



(2.1) 



provided that the following holds for all n G S n : 



trace 



H 



(n) 



> 0. 



(2.2) 



Note that in general this condition is not obviously satisfied. However if e has a 

Feynman-Kac representation, then (2.2) must be true. 

In the grand-canonical ensemble, the number of particles is not fixed, so we must be prepared 
to deal with an arbitrary number of particles. To do so we construct the Fock space J-^Ha) 
which consists of the tensor product of each canonical n-particle Hilbert space for every 
n > 0. Then the thermodynamic limit simply consists of letting the region A tend to lR d . So 



(n) 

A 



n>0 



where := C and := H A . Define the second quantization of the single-particle 
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operator A, denoted by dT(A), as the operator which acts upon J 7 (Ha) by 

dT(A) = 0#) 

n>0 

where A^ = and A^ = A. The number of particles is represented by the operator 
N A := dT(I). 

Again for the case of bosons, we need to restrict the system to the symmetric subspace of 
J 1 (V. a), denoted J-" + ("H A ), constructed by applying the symmetrisation operator to each "H A ™^ 
in the tensor product: 

F + (Ua) := n[ n l = HP = ^ E 

n>0 n>0 n>0 ' Tr£S n 

The grand-canonical partition function for a system with Hamiltonian H A defined on the 
symmetric Fock space J r + (H A ) and with chemical potential \x can be written as 

oo , 

EM = trace ^ +(Wa) [e~^^] = E^ E trace « A "> [U^ H ^] . 

?1=0 TT&Sn 

We define the corresponding probability measure on U^Lo by 

1 00 1 

(true if (2.2) is satisfied) where X is the indicator function. We take So = {1}. 

Each permutation tt & S n can be decomposed into a number of cyclic permutations of lengths 
qi, q 2 , ■ ■ ■ ,q r with r < n and gi + g 2 + • • • + q r — n. We consider the set Q = (J reN fi r of 
unordered r-tuples of natural numbers q = [qi,q 2 , ■ ■ ■ ,q r ] G ^ r for r = 0, 1, 2, ... , and let 
IqI = Qi + <?2 + • ■ ■ + q r for q G f2. Then a decomposition of 7r G 5 n into cycles is labelled by 
qGfi with |q| = n. We recall the following facts on the permutation group. 

• The decomposition into cycles leads to a partition of S n into equivalence classes of 
permutations C q , |q| = n. 
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Two permutation n' and n" belong to the same class if and only if they are conjugate 
in S n , i.e. if there exists a 7r £ S n such that 



// -i / 

71 = IX IX 71. 



(2.3) 



The number of permutations belonging to the class C q is 



V(<?i<?2 • • -q r 



(2.4) 



with n q ! = n\\ri2\ ■ ■ ■ rij\ . . . and rij is the number of cycles of length j in q. 



We observe for a symmetric Hamiltonian (i.e. [H A ,U n ] = 0,V7r £ S n ), for 7r', n" £ C q , and 
using (2.3) one has that 



trace 



n 



(n) 



U w "e 



trace ^(n) 



trace 



trace 



H 



(») 



H 



(») 



(2.5) 



Hence in the case of a symmetric Hamiltonian we may express the symmetrisation operator 
o+ as a sum (with some scaling) over all unitary representations of permutations with 
unique cycle structures of S n i.e. 



trace „,(„ 



H 



(») 



e -PH2 j = ^ ^ trace U n e~^ 

= E 



n^\q 1 q 2 . . .q r ) H ^ 



Uire 



-mi 



This motivates a more careful analysis of the class of cycle structures, specifically the ex- 
pectation of cycles of particular lengths in C q . 
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2.2 Cycle Densities 

For q e N, define the random variable N q (7r) to be the number of cycles of length q in ir. 

Fl[N q (7T)=r] = 

N q (n)=r 

Then we can define the expectation of cycles of length q in the canonical ensemble by 

n n 

K(N g ) = J2rF n A (N q = r) = ^ £ rP^[iV» =r] 

r=l r=l TreSn 

and hence the average density of particles on g-cycles for a system of n bosons is 

<5M = iSp. (2.6) 

We note here that if n/|A| = p then the sum of the densities of particles on all cycle lengths 
is exactly the density of the system, i.e. 

n 

X>M = p. ( 2 - 7 ) 

9=1 

In the grand-canonical ensemble, we similarly obtain that the expectation of cycles of length 
q is 

oo 

K(N q ) = J2rT»t(N q =r), 

r=l 

and so the average density of particles on g-cycles is 

<*) = (2-8) 

Define the grand-canonical expectation of the operator dT(A) as 

1 1 00 

(dT(A)) HA = — -.trace T+m) [dT{A)e-^-^] ^E^W^U ■ (2.9) 

n=0 
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where 

(i ( "V =^4"T trace «w U (re) e"^l. (2.10) 
-"a Z A (n) w a,+ L J 

The density of particles for finite volume is the expected number of particles per unit volume 



PaO) 



Wh a 



\A\ 

Then the density of the system is the sum of all the cycle densities: 



E^HM/i). (2.11) 
9=1 



Hence with this formulation, equations (2.7) and (2.11) imply that we can uniquely partition 
the ensemble of the particles into their respective cycle structures. 

However we would like to observe the behaviour of these cycles in the thermodynamic limit. 
To do so we shall decompose the total system density into the density of particles belonging 
to cycles of finite length (f? s hort) and to infinitely long cycles (f?ion g ) in the thermodynamic 
limit. We are thus led to the following definition: 

Definition 1 The expected density of particles on cycles of infinite length is given 
• in the canonical ensemble by 



0fong = J™ S C A^)' 

" Q— >00 n,|A|-Voo — ' 

n/\A\=p 9=Q+1 



in the grand- canonical ensemble by 



& Q— >oo A — >oo — * 
q=Q+l 
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For notational convenience, we shall henceforth define liniA as the limit as A ^ M. d , where d is 
the dimension of the system. The limited density in the grand-canonical ensemble is defined 
by />(/■*) : = lim A p A (/i). 

Note that if we define the corresponding densities of short cycles as 



f^short 



Q 

lim lim Vc^(g), 

Q— >00 n,|A|-Kx> ' ' 

n/\A\=p 9=1 



d 1 

c short 



lim lim 

)-»oo A z — ' 



9=1 



then using (2.7) and (2.11) one easily sees that 



P = Sshort + 



ong 



and p(/i) = ^ hort + q[ 



ong - 



Hence in the limit, all particles can be classified as those belonging to cycles of finite length 
and those belonging to cycles of infinite length. 



2.3 Expression for the (/-cycle Density 

Here we shall obtain a convenient expression for the density of particles on cycles of specific 
length in the following theorem. The idea is that since the trace of the partition function only 
depends upon the cycle structure of the symmetrisation operator, we can decompose this 
operator into its cycle distribution, isolate one q cycle and then recombine the remaining 
cycles into a symmetrisation operator over n — q particles leaving us with a trace over q 
distinguishable particles coupled with a boson field (with some normalisation). 

Before stating the theorem we first shall fix some notation. To indicate that a permutation 
7r contains at least one cycle of length q, we shall say q G n. 

Let U q : l-L^f 1 y be the unitary representation of a g-cycle on T-L$ defined by 

U q (4>l <S> 02 ® 03 ® • • • <8> 4>q) = 02 <B> 03 ® * ' ' ® 4>q ® 01- 
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Denote the n-space identity operator by /„, = / <g> / <g> ■ ■ • £g> I and the corresponding Fock 

n times 

space identity operator by I. When there is no ambiguity we shall simply write U q for 
U g ® I : T+CHy) -»■ jF+{H v ). Note that [U q , a\] = 0. 

Also we shall define the operator M^ A on n$ as the partial trace of e _/3( - HA ~ MArA - ) over 
JF + {H A ), i.e. 

M£ A = trace ^ +( « a) [ e ^(^-^ A )] (212) 



Theorem 2.1 TTie expectation of the density of particles on cycles of length q can be ex- 
pressed 

1. in the canonical case as 



c " (g) " zJnj\K\ trace «i 9, ®<; 9) 



(U q ® 7 n _,)e-^A 



(2.13) 



m t/ie grandcanonical ensemble as 



<{q) 



1 



trace , 



H A (//)|A| "*"~wi 9) «*+(«A) 
1 

~ 7 \ I a I trace 

=a(/x)|A| 



[(C^®I)< 



U M M 



(2.14) 



Proof: Define C n (gi,ni; q 2 ,n 2 ; • • • ; Q r ,n r ), $ distinct, q\n\ + ^2^2 + • • • + q r n r = n as the 
class of permutations with the following cycle structure: n\ cycles of length qi, n 2 cycles of 
length q 2 , etc. By (2.4) the number of elements in each class is 



#C n (q u m; q 2 ,n 2 ;...;q r ,n r ) = — 



in 



q l . . . q' r Lr n\\n 2 \ . . . n r 



From (2.1) 



K(N q = r) 



1 1 

ZJn) n\ 



trace „(») 



H 



(gi,---qk>0; qi^q\ ir£C n (q,r; q\,ny,...;q k ,n k ) 

\ Y J Qi n i= n - rc l ' 
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Using (2.5) we see that 



K(N g = r) 



1 1 

ZJn) n\ 



r; qi,n 1 ;...;q k ,n k ) 



(qi,...q k >0;qi^q\ 
V J2 qirii=n-rq ) 



x trace 



n 



-PHI 



where ir is any fixed element in C n (q, r; qi,ni, . . . ; q k , n k ). Therefore 



K(N q = r) 



#C„(q,r;qi,ny,...;q k ,n k ) 
Z A (n)n\ ' J #C„-q(q,r-l;qi,nv,...;q k ,n k ) 

fqi,-qk>o;qi^q\ 

V J2 qirti=n-rq J 



x #C n - q (q, r - 1; qi,ni] . . .;q k ,n k ) trace H („ 



1 1 



E 



??■! 



Z A (n)n\ ^ rq(n — q)\ 

V H qirii=n-rq ) 



x #C n „ q (q,r - 1; gi,rii; . . .;q k ,n k ) trace w 

^A 



We choose one g-cycle from 7r, since tt e C n (q,r; q±,n-\ \ . . . ; q k ,n k ). Thus write 7r = g o 7r', 
where 7r' G C n _ g (g, r — 1; ?i,ni; . . . ; q k ,n k ) and g G jr. The corresponding representations 
can be written as = U q ® U n > , so summing over 7r' 



PJ(JV ff =r) 



Z A {n) rq(n — q)\ 



x 



E 

V 5Z qini=n-rq I 

E 



trace 0< (n) 



ir'eC n - q (q,r~l;qi,ni;...;q k ,n k ) 



and by using (2.5) again we have 



K(N q = r) 



Z A {n) rq(n — q)\ 



x 



E 

fqi,...q k >0;qijtq\ 

V J2li n i= n - r Q ' 

E 



trace 



U 



(n) 



7r'eC n - g (q,r-l;gi,ni;...;g fc) rifc) 
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Then the canonical expectation of the number of g-cycles is found to be 



r=l 



oo 

- - y y 



X 



fqi,...qk>0;qi^q\ 

\ J2<li n i= n - r 1 ' 



trace 



H 



ir'£C n - q (q,r-l; qi,ni;...;q k ,n k ) 



1 1 



Z A (n) q(n - q)\ 



trace {U q ®U n )e~ 0H ^ 



1 1 



Z A (n) q 



trace -i/M^-i/O-?) 



(U q ® J n _,)e-^A 



and the corresponding cycle density expression follows from (2.6). Going to the grand- 
canonical ensemble, using (2.8) we have: 



<{q) 



qK(N q 



|A| EM\A\ 



n=q 



trace 



n=q 



r trace_(„) 



(/x)|A| «X"®^+C«a) 



[(17, ® I> 



-/3(i? A - M iVA)l 



as required. 



□ 



2.4 Expectation of Operators on g-cycles 



For an arbitrary one particle operator A, define the following cycle- dependent quantities: 



In the canonical ensemble, 



Z A (n)\A 



■trace , 



{A® I® 



® I)(U q ® I n _ ? )e-^A 



(2.15) 



q-terms 
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It is easily seen that c™(q, I) = c™(g). Also notice that c"(g, A) is linear in A, since the 
trace is linear. 

In the grand-canonical ensemble 



<%(q,A) = - trace (g ) 

^ A (/i)|A| H A 



(A® I® ■■■01) U q M% 



q,A 



V 

g- terms 



(2.16) 



which again is linear in A, and c£(q, I) = c^(q). 



Then the following theorem shows that we can express the expectation of those multi-particle 
operators of the form 

= A(g)I(g)---(g)I + I(g)A(g)---(g)I-\ h / ® • • • <g) / ® A 

on Tix i n the canonical ensemble, or 

dT(A) = 

n>0 

(where A^ ) := and A^ := A) on J z + {^-L A ) in the grand-canonical ensemble, in terms 
of cycle statistics. Specifically we reduce the expectation of a multi-particle operator to 
an expression where the single-particle operator is applied to one of the q distinguishable 
particles, all coupled with the boson field (with some scaling). 

Theorem 2.2 In the canonical ensemble 

M(n)\ " 

and similarity in the grand- canonical ensemble 

1 1 q=l 
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Proof: Consider the canonical case 



(A (n) ) 



— --trace vW \A {n) e~^ 
ZAn) H * 



"•A,+ 



1 1 

Z A (ri) n\ 
1 1 
ZAn) n\ 



E 

7rGS„ 

E 

TT&Sn 



trace _,( n ) 



trace 



M (n) e"^ 



using the cyclicity of the trace and the fact that [U n , if" ] = 0. Note that we can simplify 
this expression by the following method: 



trace ^ ^A {n) U^ H 

n 



trace 



trace 



i=l 
ii 



where U(n) represents the two-cycle (1 i), so again using cyclicity of the trace 



trace „,<» 



trace 



H 



(n) 



n trace 



H 



J2(A®I 

i=l 
n 

i=l 

n) \(A ®I® 

i L 



® I)U {li) U n e-^ U (li) 

® I)U {li) U n U {li) e-P H Z 
J)^e-^A 



by (2.5). Thus 

(A (n) ) 



1 1 



ZAn) (n 



-y J2 trace kW 



A L 



(A®I®---®I)U^-^ 



(2.17) 
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Define S r ({ai, . . . a r }) as the set of all permutations of {ai, . . . a r }. For some q choose 
i 2 , «3, • • • i q to be distinct elements of the set {2, 3, . . . n} and define 

S%(h,i3,---i g ) = {(1,«2, k,---i q ) °r : t e S n _ ? ({2,3, ...g}\{i 2 ,i 3 ,...ij)}. 

Using this construction, we see that the sets S^(i 2 ,i 3 , • • are disjoint for each g, since for 
each case q is the length of the cycle containing 1, the element never touched by S n _ q . Sets 
are also distinct by choice of (12,^3, ■ ■ - iq), evident when comparing two different choices of 
indices (z 2 , ■ ■ ■ iq) and (i' 2 , i' 3 , . . . i') the element 1 is in a different cycle. 

Considering the cardinality of these sets, for fixed (i 2 , ■ ■ ■ iq), #S^(z 2 , 13, ■ ■ ■ i q ) = #S n - g = 
(n — q)\ by construction. We pick the q — 1 indices (i 2 , 23, . . . i q ) from n — 1. Thus summing 
over all possibilities 

n n ( _ -\\\ 

#U U *i 3 ,..g = Ei-i# 5 n 

q=l V H >' q=l 

where S% = S%(2, 3, . . . , n). Due to the facts that S^(i 2 , z 3 , . . . z ? ) are distinct with respect to 
both the value of q and the choices of (z 2 ,«3, . . - iq), and summing all possible choices gives 
n!, we may conclude the identity 

S » = \J U ^2,...,i ? )- 

g=l ! 2 ,..i, 

For 7r G 5^(z 2 , ^3, . . . we may find cr a permutation of (2, 3, ... , n) such that 
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and therefore 



tr£ice (n) 

H A L 



(A <g) I <g) • • • <g) I)U n e 



-mi 



trace „,(«) 



trace 



(») 



trace „,(«) 



(A <g) / <g> ■ ■ ■ <g> T)U- x U£J a €-^ 
U a {A g> / ® • • • <g> 7)[/- 1 f/ # e-^A 
(^<g> J® •••<g>/)£/ # e~^A 



since U a (A <g> I <g> ■ ■ • <g> = A <g> 7 ® ■ • • <g> J and , C/ CT ] = 0. 

Thus instead of summing over all elements in S" n in equation (2.17), we may sum over all 
from q — 1, . . . , n as shown below 



1 



1)! 



trace ~, (n) 

rt A L 



(A® I ® ■ ■ ■ ® I)U^ H ^ 



1 1 - 
Jnj in - 1)! ^ 



9—1 (12,13, ■■■*«) Te5'(i2,i3v*<j) 
) |^(v4 g) / 



1 1 

Z A (n) (n - 1)! ^f(n- q)\ ^ trace «r 

9 — 1 7TI 



{A®I®---®I)U^ H ^ 



1 1 r 
^ttV/ trace („) (A <g> J <g> • • • <g> J)C/ g o U^e~^ 

l A 

2 , trace ^ 

9=1 



(A <g) / <g) • • • <g) I)U q e~ 



as desired. The grand-canonical case follows with ease. 



□ 



In the following chapters we consider Feynman's conjecture for a selection of models and 
investigate the relationship, if any, between a non-zero infinite cycle density Qi ong and Bose- 
Einstein condensation. 
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Cycle Statistics of Simpler Models 



Summary 

In this chapter a relation between the ground state occupation density and cycle statistics is 
made. This relation rests upon the validity of three conditions, all of which are verified for 
three different models: the ideal Bose gas, the Mean-Field model and the Perturbed Mean- 
Field model. From this theory we reproduce the well known fact that the density of long cycles 
for these three models is equal to the density of the condensate. However we note that the 
core condition does not hold for all models. 



3.1 Occupation Numbers 

The appearance of (generalised) Bose-Einstein condensation depends on the low-energy be- 
haviour of the density of states. We intend to investigate the connection between the cycle 
statistics and the occurrence of condensation by considering the expansion of the expectation 
of a low-energy density operator in terms of the formulation described in Section 2.4, and 
comparing it with the long-cycle density £>i on g- 

Consider a system of bosons enclosed a region A c R d . Let 

h A = -A 
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act upon a dense domain in the Hilbert space l-i A : = L 2 (A) with suitable boundary conditions 
on dA. Suppose that h A has eigenvalues Ef; with corresponding eigenvectors 0£, indexed by 
k > 1 such that < £ A < £ A < . . . . 

Define a symmetric Hamiltonian H A on the Fock space J-+(H A ) by 

H A = dT(h A ) + U 

where U denotes the inter-particle interaction. For any state G % A define the orthogonal 
projection in H A onto <f) by P^ := The corresponding second quantization is the 

operator N A (4>) := dT^P^) which counts the number of particles in the state 0. Let P £ := 
Y2k-E A <e P<t>k w hich is a projection onto the subspace of % A with single particle kinetic energy 
not exceeding e. Then the operator 

N A {e) = dT{P £ )= N M) 

k:E£<e 

counts the number of particles with kinetic energy less than e. 

In terms of cycle statistics established in Chapter 2, the average expectation of this operator 
can be expressed by 



(N A (e)) HA 



" v ' 

g-terms 



where is a partial trace of e~^ ,H ' A ~ AtJV ' A ) over J 7 + ('H A ), see (2.12). We shall consider this 
cycle expectation expression to investigate whether cycles of infinite length appear and if 
they coincide with a macroscopic occupation of particles near the ground state. 

Using this formulation, under certain conditions, we can prove three Theorems relating the 
density of the Bose-Einstein condensate with the long-cycle density, which can be applied to 
the simpler models discussed in the next section. These Theorems depend upon the three 
conditions, which we shall number I — III. Later we shall prove the validity of these conditions 
for a selection of models. 
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For the next three Theorems let us first assume that the following condition holds: 

(%(q,A) > for all A > 0. [ I. ] 

The first theorem proves that the occupation density for any single particle state is bounded 
above by ^ ong . 

Theorem 3.1 Suppose that Condition I is satisfied and if 

limc£(g,P, A ) = (3.1) 

holds for all q > 1 and any (fr A G H A with ||0 A || = 1. Then for any single-particle state 
4> A £ 'Ha we have the following upper bound: 



W )) 

|A| 



Hp Al^< Cg . (3.2) 



Proof: 



oo Q 

lim <^ I~P<P A )> lim lim £ c^g, I - P 0A ) 



lim lim - lim lim V c£(g, P 0A ) 

Q— >oo A z — ' Q-s.oo A 

9=1 9=1 



lim limVc^(g) 

O— ior, A ^— » 



¥oo A 

9=1 



using Condition I and (3.1). Then since lim^^c^(g) = p(/i) and c A (q, A) is linear in A, we 
may write 



A 

9=1 



oo oo 



lim^^^lim^c^P.J^^lim £ <(<?) = < ong 

I I 9=1 °° 9=Q+1 

as desired. □ 
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We note that condition (3.1) holds if we assume that e~@ H has a Feynman-Kac representation. 
In this case we can show that 

c£(?,P,J < < (3.3) 

where E$ is the lowest eigenvalue of the single-particle free Hamiltonian h A , which immedi- 
ately implies condition (3.1). 

Justification of (3.3) is as follows: let be the restriction of H A to Ji A n \ The kernel of 
exp(— /3iJ A ) on 7i A ^ ®l-L A n ^ is positive and bounded above by the kernel of exp(—/3Hl) (8 
exp(—/3H2~ 9 )- Thus the kernel of M^ A is positive and bounded above by S A (/i)e^ M x the 
kernel of exp(—/3/i A ^). Therefore 



1 



<(g, < j^jtrace^) [(P lM ® / ® . . . ® t)U^T^-^ 
Then considering the trace and expanding in terms of the eigenbasis of /i A 



trace ^(a) [(P|0 A | <8 I® ■ ■ ■ <8> /)!7 g e" 



E e" /3( ^i + ^ + - + ^ ) 

fcl,fe 2 ,...,fcg>l 

X <0fe x ® < ® • • • ® <f> A k q \(P\<p A \ ® / ® • • • ® ® 0fc 2 ® • • • ® O 

^ e -^(<+^ a +-+^) 

fcl,fc 2 ,...,fcq>l 

x (0* ® 0£ 2 g) . . . ® ^J(P| 0A , (8) / (8) . . . (8) /)|0t 2 ® 0fc 3 ® • • • ® O 

fci,fc2v;k g >l 

= ^e-^ A <^||0 A |)<|0 A ||^> 
fc>i 

< Q -PqE 1A 



and hence (3.3) follows. 

In the next theorem we show that the generalised condensate in the lower single particle 
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kinetic energy states, defined by 



r ,. (N A (e)) 



e^O A |A| 

is bounded above by Qy ong . 

Theorem 3.2 // Condition I is satisfied and 



lim lim c^(g,P £ ) = 0, [II. 



then for any single-particle state A G H A we have the following upper bound: 

^cond — ^long 



Proof: Similarly to the previous proof we have that 



lim V c£(g, I-P e )> lim lim V c£(q, I - P e ] 

A z — ' Q->oo A z — ' 

9=1 g=l 



Q <3 
lim lim V<((|) - lim limVV(g,P £ ) 

Q— >oo A ^— — ' Q— 5>oo A z — ' 
g=l g=l 



lim limVc^(g) 



>oo A 

9=1 



using Conditions I and II respectively, with which one can derive the following: 



oo oo 



lim lim V. . = lim Y <{<li p e) < lim lim Y c a(<?) = C- 

1 1 9=1 g=Q+l 



Finally we show that is also an lower bound for the generalised condensate. 



□ 
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Theorem 3.3 Suppose that for each e > 0, 



then 



lim lim V d&q,I-P e ) = 0, [III. 

Q->oo A * — ' 
q=Q+l 



^cond — ^long 



Proof: Let Q £ = I - P £ and recall that (N A (e)) := (dT(P B )). Then 

oo Q oo 

lim V c»{q, Q e ) = lim lim V c£(g, Q e ) + lim lim V c^(g,Q £ ) 

A z — ' Q-s>oo A z — ' Q->oo A * — ' 

9=1 9=1 9=Q+1 

Q 

= lim limVV(g,Q e ) 

Q— >oo A z — ' 

9=1 

Q 

< lim limVc^(g) 

Q->oo A z — ' 
9=1 

by Conditions I and III. Hence 

limlim — ^— — = lim lim > eft (a, P e ) > lim lim > dt(q) = o{^ no . 

e^O A A e^O A ^ AV ' ~ Q^oo A ^ AV ' long 

1 1 9=1 9=<2+l 



□ 



Consequently those models which satisfy Conditions I — III, i.e. the foundations of Theorems 
3.2 and 3.3, are endowed with the property that 

ftong ^cond' 

We shall prove these conditions for three models, the Ideal Bose Gas, the Mean-Field model 
and the Perturbed Mean-Field model. In the following section we give a brief introduction 
to these three models, and prove the conditions in Section 3.3. 
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3.2 The Models 



The regions for these models, A;, are a sequence of d- dimensional cubic boxes centred at 
the origin with sides of length / = 1,2,..., i.e. A; := [— |, and thus with volume 
Vi : = | A; | = L d . We shall relabel all quantities with / instead of A. Thus the single-particle 
Hamiltonian h\ = —A acts upon %i := L 2 (h{) with Dirichlet boundary conditions on the 
boundary, and we label its eigenvalues by < E\ < E l 2 < . . . with corresponding eigenvectors 

The ideal Bose gas Hamiltonian is the second quantization of hi on the bosonic Fock space 

= dT(hi). (3.4) 

The thermodynamic limit is formed by letting / — > oo, i.e. Vi — > oo. The phenomenon of 
Bose-Einstein condensation occurs when eigenstates of the Hamiltonian are macroscopically 
occupied in the thermodynamic limit. To examine this behaviour, we define a distribution 
function vi which counts the number of eigenvalues per unit volume by 

*i(A) = ^#{k : E[ < A} 

vi 

for any A G K + . One can prove (see [93] and [89] for more details) that this sequence 
of measures converges weakly in the thermodynamic limit to a measure u, known as the 
integrated density of states, which in this case evaluates to 

u(X) = C d \ d/2 (3.5) 

where C d = (vr d / 2 2 d r(f 

In terms of the single-particle eigenstates <pk (dropping the / for clarity), define the usual 
creation and annihilation operators a*(0fc) and a{(j>k) which satisfy the canonical commut- 
ation relations [a(<j>j), a*(<pk)] = $j,k- The operator Ni((j>k) '■= a*(4>k)a{.4>k) then counts the 
number of particles in with state 0^. Using this formalism we may rewrite the 
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grand-canonical Hamiltonian of the ideal Bose gas, equation (3.4), as 

fl?(Ai) := Hf - /iiV, = J2( E l ~ A*) Nitek) 

k>l 

where /i is the chemical potential of the ensemble and iVj := ^ fe>1 -Nz counts the total 
number of particles. 

For the ideal Bose gas in a finite volume with chemical potential \i < 0, at thermal equilibrium 
with inverse temperature j3, we define the grand-canonical partition function by 

S°(/*) = trace ^ +fMi) [ e "/»W] . (3.6) 

This allows us to define the finite- volume pressure and density of the ensemble as 

ft°(A0 = ^jh3°(A0, rf(A*) = 

We denote the corresponding limiting quantities by := limj_ ) . 00 £>j ) (jt/) and p°(yu) := 

lim/^oo p°(/x) respectively. 

The expectation (Ni(<fik))H°(ij,) returns the average number of particles occupying a specific 
quantum state <pk when the system is in thermal equilibrium. It may be evaluated to yield 

and hence we can re-express the total average density of the system in finite volume as 

= y t £< W fc )>H» = / ^TTT^(A). (3.7) 

Similarly one can derive the following expression for the pressure 

ftV) = 4 / ln[l-e-^]^^(A). (3.8) 
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It is useful to define the partial pressure n : (—00, 0) —> R + by the following 

7 r( S ) = iln[l-e^]- 1 (3-9) 

which allows us to write 

p?( A *)= / ttGu-A)^(A). 

To observe the well-known phase transition is this model, we proceed as follows: fix the 
mean-density of the ensemble as p. We then require that the chemical potential p satisfies 

pKw) = p 

for all I. Set p^ : = lim/^oo pi, and denote the critical density by 

Per = lim p{p) 

and hence 

Pcr = [ (e^-l)- x dKA)- (3.10) 

Note that the limiting thermodynamic quantities like p(p) only exist for negative values of 
p. But for more general models, p can be positive or negative depending on the interaction 
and the boundary conditions of the system. 

If p cr is infinite, then we find that p^ < always satisfies the equation 

P = p\poo) (3.11) 

which implies a lack of a phase transition (this occurs for d = 1,2, see [49]). If p cr is 
finite (here d > 3), then equation (3.11) is again solved by a unique non-positive p^ when 
p < p cr . However if p > p cr then (3.11) has no solution and p^ = 0. It can be shown 
that the corresponding grand-canonical pressure p°(p) is constant when p > p cr , indicating 
a phase-transition at p cr to which we associate with the appearance of generalised Bose- 
Einstein condensate. To locate this condensate we focus on the eigenstates with energies 
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in the neighbourhood of the ground-state before taking the thermodynamic limit. For any 
e > let 

N t (e) = WW 

which counts the number of states with energies less than e. Then one may calculate the 
average number of particles with arbitrarily low energies to find g COQ d, the density of the 
condensate: 

&ond := hm hm = Inn hm / -dvi{\) 

r i z" 00 i 

= p — lim lim / —r- — r -dui(X) = P— — 577 ; -dv{\) 

H e^oi^oo J [et0o) e^( A "«) - 1 Jo e^( A -^-)-l v ; 



f?cond 



P < Per 

P - Per P > Per 



(Here there is no p dependence as we fixed the mean density of the system p. In the other 
models we consider this will not be the case so we shall denote the condensate density by 
PcondO ^ ^ s fascinating that a collection of non-interacting bosons can yet display a phase 
transition, purely as a consequence of their quantum mechanical properties. As a result, 
this model has been the object of study for many years, and has been approached in various 
ways, for examples see [52; 72; 90; 92; 94]. 

Note also that here we did not consider the nature of the condensate, but one in fact can 
show that the ground state is macroscopically occupied 

MOM) - 

lim 77 = P ~ Per- 

I— too Vi 

However it has been shown by van den Berg, Lewis, Pule [93] (see also [86; 90; 92]) that the 
ground state is not necessarily the only state to be macroscopically occupied, but instead 
observed there are three possible ways the condensate may form, depending on the shape of 
the region A;. Instead of taking A; to be a d- dimensional cube, they considered rectangular 
parallelepipeds whose edges go to infinity at different rates (known as Casimir boxes [21]). 
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A consequence of this anisotropy is that the standard ground-state BEC is converted into 
generalised BEC of three types. These condensates may be classified as follows: type I has 
a finite number of eigenstates macroscopically occupied, type II an infinite number of states 
occupied and most surprisingly a type III, which has no macroscopically occupied eigenstates 
before the thermodynamic limit, but where a macroscopic number accumulate in the ground 
state in the limit. 

This thesis will only deal with the generalised condensate, but an investigation of the relation 
between cycles and these three condensate types is currently being performed by M. Beau 
[6]. 

With respect to cycle statistics, Siito showed in his two papers on the topic [78; 79] that 
£cond exactly equals the long cycle density g£ for the ideal Bose Gas. 



The Mean-Field Model 



In the mean-field model of a system of interacting bosons (also known as the Imperfect Bose 
Gas), an external energy proportional to the term Nf jV\ is added to the Hamiltonian of the 
ideal boson gas: 

H f = H? + —Nf 
i 1 2V t 1 

where a > controls the strength of the interaction [50]. 

This model has been studied extensively, see for instance [24; 87; 71; 91], where it is shown 
that Bose-Einstein condensation persists in the presence of this external interaction. Con- 
densation only occurs if it occurs in the case a = 0, i.e. the ideal case, with similar depend- 
ence on the dimensionality. Note that the presence of the mean-field interaction stabilises 
the system so that the thermodynamic quantities are defined for any p. Here the condensate 
density is 

{0 p< ap cr , 
p cr p> ap cr . 
a 
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where p cr is the critical density of the ideal Bose gas. 

However in contast to the ideal Bose gas, the occurrence of the condensate is due to an 
effective attraction between bosons in the zero-mode [18]. Moreover at fixed mean-density 
it is independent of the strength of the repulsive interaction and equal to the total amount 
of condensate in the Ideal Boson gas at the same mean-density. 

Siito again proved the equivalence of long cycles and the condensate for this model. 
The Perturbed Mean-Field Model 

This model is an extension of the mean-field model with the addition of some interactions 
between modes. Its Hamiltonian is 



where Nk '■= Ni(4>k) is the occupation number operator of the mode k, a > and we write 
v(k,k') to denote v(E] c , E^r), where v has a continuous and positive definite kernel with 
suitable decay properties. 

Originally considered in 1990 by van den Berg, Dorlas, Lewis and Pule [87] (c.f also [33]), 
they prove the existence of the thermodynamic pressure and find a variational formula for 
this pressure using Large Deviations (we shall lay out the rough procedure of this paper as 
the method will be required later). Later in 2005, Dorlas, Martin and Pule [34] consider the 
cycle statistics of this model and prove that the density of long cycles equals the condensate 
density. 

To apply Large Deviation Theory, one must first phrase the question in terms of probabilities. 
The Hamiltonian (3.12) is diagonal with respect to the occupation number operators Nk, so 
it is possible to consider the occupation numbers as random variables rather than operators. 
Define our probability space Q as the countable set of terminating sequences of non-negative 




(3.12) 
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integers, so that any u £ f2 is of the form {u\,U2, ■ ■ ■ \ oo j £ N V? with^Wj < oo}. Then 
the basic random variables are the occupation numbers {iV 1; N 2 , . . . }, which are maps Nj : 
Q — > N defined by Nj(u) = Uj for each configuration u; £ f2. The total number operator 
iV(u;) for a configuration a; is defined by 

i>i 

which is necessarily finite by construction. With respect to this formulation, the Ideal Bose 
gas Hamiltonian in the region A/ with eigenvalues < E[ < E\ < . . . (henceforth neglecting 
the I) may be expressed for a configuration u as 

H°(u)=^E j N s (u) 

and thus the Perturbed Mean-Field Hamiltonian is 

aN{uf+ Yl v(k,k')N k (u)N k r(u) . 

k,k'>l 



Hr f (u;) = H?{u) + — 



In this representation the corresponding finite- volume grand-canonical partition function and 
pressure for a chemical potential \i £ R are (respectively) 

Again we set the limited pressure p(/i) = lim^oop^/i). Recall for \i < that equations (3.6) 
and (3.8) define the ideal-gas partition function H°(/i) and pressure p°(/i). 

Using Large Deviation theory, a variational expression for p(/j) can be found with respect to 
the pressure of the ideal Bose gas. We begin by placing a probability measure upon Q with 
respect to H®: 
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Introduce the occupation measure L\ by defining for each Borel subset A of [0, oo) and u G Vt 

Li[u;A] = y^N j (u)5 Ej (A), 

where S x (-) is the Dirac measure concentrated at x. Then for every u in Q, L\ is a bounded 
positive measure. It is convenient to define the following notation for any m G .Mi! (R+), 
the positive bounded measures on [0, oo): 

\\m\\= I m(d\), ((m,Vm))=[[ v(X, X')m(dX)m(dX'), 

and put 

G M [m] = n\\m\ | — — | \m\ | 2 — — ((m, Vm)) . 

In terms of the probability measure we can express the pressure for the interacting 
model as a ideal gas pressure plus a correction term which incorporates the inter-particle 
interactions (given by G M ): 

W (/i)=rf(«) + ^ 7 ln^e^- Q ^[--¥r(c) (3.13) 

for a < (since the ideal-gas pressure is only defined for strictly-negative chemical potentials, 
we use a trick from [89] whereby negative a is used for this purpose, for which we compensate 
by introducing a \x — a term in the interaction term G^ a ). 

Next we rewrite (3.13) as an integral over £ := Ai b + (M. + ), which we equip with the narrow 
topology to guarantee the continuity of the mapping m — > G M [m] (see Lemma 4.1 of [87]). 
G M also has an upper bound since v is positive definite, ((m, vm)) is non- negative and a > 0. 
Let K/ be the probability measure induced on £ by Lf. 

K? = P?oLf\ 
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with which (3.13) can be re-expressed as 




(3.14) 



The form of the above integral suggests the use of Laplace's method to evaluate it. In fact, 
we shall use a variant called Varadhan's theorem [95], which depends upon the measures Kj 
satisfying certain criteria so that for large Vi, they exhibit the following behaviour 

K/[dm] ~ exp ( — Vil[m\)dm 

where I[m] is called the rate function. The criteria are stated as follows: we say that the 
sequence of measures K/ satisfy the Large Deviation Principle [38; 95] with constants V\ 
(which diverge to +oo) and rate function I : £ — > [0, oo] if the following hold: 

(LD1) I is lower semi-continuous, 

(LD2) for each b < oo, the set {m : I[m] < b} is compact, 
(LD3) For each closed set C 

lim sup — In [C] < — inf{/[m] : m G C}, 

l^oo Vl 

(LD4) For each open set O 

liminf — InKjO] > - inf{J[m] : m G O}. 

I— >oo Vi 

Now we may state Varadhan's Theorem: 
Theorem: Varadhan [95] 

Let Ki be a sequence of probability measures on the Borel subsets of £ satisfying the Large 
Deviation Principle with constants V and rate-function I. Then, for any continuous function 
G : E — > K which is bounded above, we have 

lim l~ In ( e^ G[m, Ii[H = sup {G[m\ - I[m]} . 
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It is shown in a paper by R.S. Ellis, J. Gough and J.V. Pule[39] that there is a large class 
of measures with random weights which satisfy the large deviation principle and they give a 
general procedure to obtain the corresponding rate function in an explicit form. It is easy to 
check that IK" is of this class, so by taking the thermodynamic limit, Varadhan's Theorem 
implies that (3.14) reduces to 

p(p) = p°( a ) + sup [G a ~^[m] - I a [m]} . (3.15) 

Then to find an expression for the rate function J a [m], Ellis, Gough and Pule describe the 
following procedure: first single out the part of the measure m which is singular with respect 
to the integrated density of states du, i.e. let m = m s + m a be the Lebesgue decomposition 
of m with respect to dv, where m s is the singular part and m a the absolutely continuous 
part. Let p be the density of m a , i.e. m a (dX) = p(X)du(X). Define U : £ — > [0, oo] by 

U a [m] = - [ (a- X)m s {dX). 

J[0,oo) 

Let ir* : R — > (—00,00] be the Legendre-Fenchel transform of ir, that is 

7r*(t) = sup {ts — vr(s)} 

where tt is defined in equation (3.9). 
For t G IR and r < 0, denote 

J(t,r) = 7r*(t) -rt + n(r). 

Then Theorem 3 of [39] gives an explicit expression for the rate function of the form 

r[m] = U a [m s ] + f J(p(A), a - X)m a {dX) 

= - [ (a - X)m s (dX) + I Tt*{p{X))dv{X) + p°(a) 

J [0,oo) J[0,oo) 

and hence from (3.15) get an expression for the limiting pressure in terms of a variational 
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problem 

= - inf {(Ei-/i)||m|| + ^||m|| 2 + ^((m,1/m))+ / ir*(p(\))dv(\) 

rn&£ I 2 2 J[0,oo) 

The existence of a minimizer is proved by Dorlas, Lewis and Pule [87], but an analytical 
expression remains to be found. However they succeeded in showing that if v (•, •) is continu- 
ously differentiate and dv/dX is bounded, then the minimising measure has an atom at the 
origin, i.e. m s is concentrated at A = 0. With this Dorlas, Martin and Pule [34] prove that 
the density of long cycles equals the condensate density. We shall reproduce this result in 
the next section by a different method. 



3.3 Checking the Three Conditions for Simple Models 

The technique detailed in Chapter 2 depends upon finding the probabilities of particles on 
cycles of particular lengths. The validity of these probabilities rests upon (2.2), ensuring for 
the model chosen that the Hamiltonian Hi satisfies 

trace (n) [U n e~ m ] > 

for any n£l, and for all ir £ S n . 

We shall begin by checking this for the Perturbed Mean-Field model. As will often be 
the case, in considering the Perturbed Mean-Field model, the result will also stand for the 
Mean-Field by taking v — 0, and the Ideal Gas with a — v — 0. 

To prove (2.2) for the Perturbed Mean-Field model, we first shall denote the basis of unsym- 
metrised Hilbert space nl n) by |k) = |A;i, . . . , k n ), for all k £ N n . Then it is easy to see that 
N\k) = n|k) and A^|k) = rij(k)|k), where nj(k) £ {0,1,..., n} such that Ylj^^jO*-) = ti. 
We expand the (unsymmetrised) trace in terms of this basis as follows: 



trace _(„) 



U n e~^ = k 2 , . . . , k n \e^ H i m U v \h, k 2 ,..., k n ) 
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k 

{71 1 " "I 



|^7r(l)j ^7r(2) 5 • • • j ^7r(n)) 



exp <^ -/3 



2V^ 2^ 



1 n 

Y7 X/ u (i>/K( k K'( k ) 



which is a sum of exponentials of real numbers, the resulting sum of which must be non- 
negative. 

With this foundation in place, Theorems 3.1 and 3.3 state that if the symmetric grand- 
canonical Hamiltonian Hi satisfies the following conditions: 

I. cf (q, A) > for all A > and q e N, 
II. for any gel, lim lim cf(q, P e ) = 0, 

e— >0 I— 5>oo 



III. for some £ > 0, hm hm > c f (?> ^ — -Pe) = 0, 

Q— >oo I— >oo — ' 
g=Q+l 



then the model is endowed with the property that the density of cycles of infinite length 
equals that of the Bose-Einstein condensate, i.e. 



^lone ^cond 



We shall consider the following three models introduced above with this framework: the 
Ideal Bose Gas, the Mean-Field Model and the Perturbed Mean-Field Model, to reproduce 
the results of Siito [78; 79] and Dorlas, Martin and Pule [34], i.e. that the condensate density 
equals the long cycle density. Again it suffices to assert Conditions I — III for the Perturbed 
Mean-Field Model, since setting v — returns the Mean-Field Model, which also implies the 
result for the Ideal Bose Gas. 
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However before considering these three conditions, we shall first state and prove a useful 
lemma. 

Lemma 3.1 For any non-negative single-particle operator A we have the following identity: 

■ """V k>l ly k^<J 

Proof: Before we begin we need to establish some notation. Denote an ordered set of q (not 
necessarily distinct) positive numbers by k := [k\, k 2 , . . . , k q ) G N q and let 

|k) = 0fci ® <Pk 2 ® ■ ■ ■ ® <Pk q - 

Then {|k)|k G N 9 } is an orthonormal basis of W^f 1 , and a summation over k will sum over 
all possible g-tuples. 

Let {ip n }^Li be an orthonormal basis for J z + (1-Li). Then a basis for T-t^ cg> F + {Hi) may 
be formed by taking the tensor product of the bases of %i and J z + (l-ii), thus the set 
{ |k) ® ip n := |k, ip n ) | k G W; n — 1, 2, . . . } is an orthonormal basis of <8> ^(Hi). 

Since the expression for the grand-canonical cycle density, equation (2.14), involves a Hilbert 
space of the form nf^ <S> J-*+(K/), it is convenient to establish a notation for operators acting 
on T-L^ and J r + (jH,i) individually. 

Define N% on as the operator which counts the number of q unsymmetrised particles in 
the state 0^. This may be expressed as 

q g 
Nf\k) = N?\(ki,k 2 ,..., k q )) = kkMh, h,..., k q )) = |k> kki = |k)^(k). 

3=1 3=1 

where n»(k) G N. Note that Yli>i n i(k) = Q- The q- unsymmetrised free-particle Hamiltonian 
is hi which has eigenvectors |k) for all k G N 9 , with corresponding eigenvalues 

i=l 



cf(q,A) = — — -^2(0 k \A(p k ) trace r+iHt) 
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Therefore the operator corresponding to N k applied to 1-if' <8> F+il-Li) is + N k , and the 
total number operator is q + N. 

For the remainder of this proof we discard the I subscripts from the Hamiltonian H, the 
volume V and the number operator N for clarity. By definition 



c?(q,A) 



1 



-trace 



(A <g> I <g> ■ ■ ■ ® I)(U q <8> l)e-^ HPmf -^ 



which when expanding the trace in terms of the basis of Tli <8> jF+iJii) becomes 



= T/-f ^ ^(k,^ w |(A®/®---®/)(t/ g (g>]I) 

k 71=1 

{1 11 
-0 GT(/l) + ^iV 2 + — J2 V ( k ^ k j) N i N J-f^ N \\ k ^n)- 
i,j=l -I J 

Now split the operators Nj and H° into their "Hp ® J r + ('Hi) representations, Aj + Aj and 
/i^ + dT{h) respectively, to obtain 



1 °° f 

^— ^ ^ <k, V>„ | (A (8) / <g> ■ ■ ■ <g> I)(U q ® I) exp <^ - p 



mil*) 



k n=l 



^) + dr(/i) + _ {q + N) 2 



OO -, s 

+ 27 E "(^J) W + Ni)(N] + iV.) _ fa + AT) L | k , ^} 

i,j=l - 1 J 



All the operators in the exponential commute, and as A^|k) = n/t(k)|k), this is equal to 



k n=l 



■•® /)(£/, ® I) expj -/3 
^ ]T j)(n 4 (k) + NMnjQt) + iV,) - + A) 



^) + dr (/ i ) + — (g + iV) 2 



|k,^ n ). 



Since 



(k|L4®/® ••■(g)/)C/ 5 |k) = 



unless k\ = k 2 = . . . = k q := k, we can write the same expression over one k and replace the 
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rij(k) with q5i. as follows: 



< (g, A) = "£{k\Ak) E W ex P I - /3 



n=l 



^) + dT(/ i ) + — (g + iV) 2 



j oo 

— ^ v(i,j)(q5 ki + Ni)(q6 kj + Nj) - fi{q + N) 



*>J=1 



|^n) (3.16) 



- — . . 'S^lklAk) trace 



N k >q 



exp < — /3 



l j>i 



+ —N 2 , 
2V 2V 



j oo 



— — E(^I AA; > trace 
fc>i 



/3{H<> mf -ij,N) 



where we mean that the trace is taken over the subspace of Fock space with at least q 
particles in the state (j) k . D 

Corollary 3.1 Let <p k be an eigenstate of the ideal Bose gas Hamiltonian, and P^ k the 
orthogonal 'projection upon this state, then 



1 



trace 



ViZi(fj) Nk > q 



e -^(flf mf _ M JV,)' 



Now we may proceed to consider the three conditions proposed in Section 3.1. 



3.3.1 Condition I 



Now we shall check the non- negativity of cf(q, A) for any single-particle operator A > for 
the case of the three models described above: the Ideal Bose Gas, the Mean Field model and 
the Perturbed Mean-Field model. By Lemma 3.1 we have that 

e -/3(fff mf -^)" _ 

k>i "k^q 



cf(q,A) = ^{(j) k \A(j) k ) trace T+m 
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Since A > 0, then (<f)k\A(f)k) > for all k > 1. The trace term in the summation is positive as 
it involves the exponential of positive operators. Hence df(q, A) > for any q and Condition I 
is satisfied. 



3.3.2 Condition II 



Here we wish to prove for the Perturbed Mean-Field model that 



lim lim cf(q, P e ) = 



where P £ := 'Ek:E k <e P w 



By Lemma 3.1, for a state (f) e ~H A we immediately obtain 



tf(q,P<i>) = , r „, x Y] 1 (0 1 0fc) 1 2 trace jr m 



-/3(//, pmf -^) 



The trace with the restriction Nk > q is bounded above by the same trace but with > 0, 
which actually is the usual partition function for the Perturbed Mean-Field model. Hence 



1 



trace F+{Hl) 



N k >0 



1 k>l 



which tends to zero as / — > oo. Then using this 



lim lim cf(q, P £ 



lim lim cf(o, P^, ) < lim lim — 1 



k:E k <e 



e— >0 i— >oo V; 



k:E k <e 



e— >0 I— >oo 



< lim lim / dui(X) < lim / <iz/(A) 



[0,6) 



£->0 



[o, £ ) 



2C d 



[0,6) 



d + 2 e->o 



lim^ 2+1 = 0. 



for any d G N, recalling (3.5), the expression for the limiting integrated density of states. 
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3.3.3 Condition III 



Here we shall prove that for some e > 0, with P £ := Ylk-E < e -^W that 

oo 

lim lim > c f (OS I ~ Pe) = 0. 

Q— s>oo £— »oo ^— ' 
q=Q+l 

It is convenient to note that we may interchange the summations using Fubini's Theorem, 
i.e. 

oo oo 

£ 4(q t I-P<)= £ £ cf(g,P,J (3.17) 

<?=Q+1 fc:B fe >eg=Q+l 

since the total sum is bounded above by the total density p([i). 

The approach for the Perturbed Mean-Field case will involve Large Deviation theory. Cer- 
tainly while this calculation will imply that Condition III holds for the case of the Mean- 
Field model too, there is a short and intuitive alternative proof for the Mean-Field case alone 
which we shall give. Also to address possible concerns with guaranteeing the negativity of 
the chemical potential for the Ideal Bose Gas, we will consider each model individually. 



Ideal Bose Gas case 

Due to the absence of interactions between particles, it is trivial to deduce from equation 
(3.16) the following expression for the cycle expectation of an operator A > 0: 



for \i < 0, and hence 



Vi 

1 k>\ 



Then taking (3.17), one immediately obtains 

°° 1 OO -. MQ +1 \ E 

E E^-^-v, £ £ 4 E Wr 

q=Q+lE k >e k:E k >e q=Q+l k:E k >e 
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For large enough I there exists an e > satisfying E^> e such that 



E 



y ^ e ?E k _ i 

1 k:E k >e 



-/3(Q+l)E k e -p(Q+i) £ 
< 



E 



y e l3E k _ i 

1 k:E k >e 



[e,oo) 



e /3A _ ! 



n(dA). 



In taking the thermodynamic limit, one thus obtains the upper bound 



oo „ 

E E ^ e_W£ / ^-T-( rfA ) ^ ^ 



which tends to zero in the limit Q oo. 



Mean-Field case 



One can easily deduce an equivalent form of Lemma (3.1): 



-4- trace , +m [e^-M">]. 



Expand the Hamiltonian and incorporate the N^> q restriction to obtain 



cf(q,P^ k ) = trace F+m 



exp J-/3 J2( E i + S 3*Q) N 3 + + g)2 " (iV + 



e -Pq{E k -E!) 



-trace 



*+(«») 



exp < — /3 



+ ^7(iV + g) 2 -(iV + g)/i 



e -^ ? (S fc -Si) 



trace [e 



-/3(if" mf -^jV)] 



< 



-PqiEk-E!) 



trace [e^"^] 

iVi>0 



e -Pq(E k -Ei) 



Then taking (3.17), for large enough I there exists an e > such that Ei < e/2 and 



oo 

E E <?(^J< E 



k:E k >eq=Q+l k:E k >e q=Q+l 



-l3qE k 



E 



-PQE k 



V ^ e PE k _ I 
k:E k >e 
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which in taking the thermodynamic limit, is bounded above by 



-/3QE k 



lim 

i^oo y L^t e ftE k _ i 

k:E k >e 



which goes to zero in the limit Q — > oo. 



Perturbed Mean-Field case 



Recall Lemma 3.1 gave us, 



cf(q,P 



l 1 <t>kj 



trace 



ViEi(ti) Nk > q 



-/3(Hf mf -ixN,) 



For all fieN, define the projection P(N k = n) which projects onto the subspace of J : + (1-Li j 
with exactly n particles in the k th mode. With this we can write 



E 



n=q 



trace 



Then performing the q summation of (3.17) first we have 



oo oo 



trace P(N k = n)e 



q=Q+l n=q 



and by rearranging the sums, get 



1 

i 



oo n 

trace 



n=Q+l g=Q+l 

oo 



P(N k = n)e 



^ (n - Q) trace P(N k = n)e 



n=Q+l 



^2(n - Q)Q{n - Q - 1) trace jr +(Wi) P(N k = n)e 



n=0 



- ^ <(n - Q)6(n - Q - l)P(jV fc = n)) Hr 



n=0 
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-((N k - Q)Q(N k - Q - 1)) 



V, 



H 



pmf 



where O is the Heaviside step function. With this (3.17) becomes 



cf(g,/-P £ ) = - ((N k - QMNk ~ Q ~ 1)) Hr < 



q=Q+l 



k:E k >e 



i ®(E k - e)((N k - Q)Q(N k -Q- l)) Hfmf 



k>l 



\ \J2( N " - Q) Q ( N ^ -Q- i) ^ - £ ) 



,k>l 



H 



p m f 



This inspires us to define the following modified "free" partition function 



trace 



exp { -p H° -fiNt-r 



J2(N k - Q)Q(N k -Q- l)Q(E k - e) 



k>l 



and denote its corresponding pressure by pKfi, r, Q, e). Similarly, define a modified Perturbed 
Mean-Field partition function 



trace T+m 



exp { -P | H? mf -fiNi-T 



J2(N k - Q)Q(N k -Q- l)e(E k - e) 



k>l 



with pressure pi(fi, t, Q, e). These pressures must satisfy the following relation: 



d_ 



Pl (li, t,Q,e) = i ( V(iV fc - Q)Q(N k -Q- l)Q(E k - e) 

r=0 V l 



,fc>l 



:= <t(Q>I-Pe)- 

jjpmf q=Q+l 



We shall apply the technique of Large Deviation Theory to write pi(n,T,Q,e) in terms of 
the "free" pressure p° r, Q, e) and find its asymptotic behaviour, which will allow us to 
derive an expression for the derivative of pi(fi, r, Q, e) wrt r. 
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The "free" partition function can be written as 

1 i _ e -P(E k -»-T) + _ ^e-KQ+lX^-/,) 

k:E k <e k:E k >e V /V ' 

We can rewrite the limited "free" pressure as 

p°(/x = a,T,Q,e) = / 7r(a; - X)du t (X) + / 7r riQ (a - X)dui(X) 

J[0,e) J[e,oo) 

with r > 0, where 

1 _ e /3(s+r) _|_ ^/3r _ -^ e /3s(Q+l) 



^(s) = (3 In 



1 -e^)(l -e« s +^)) 



(3.18) 



for s < min{0, — r}. Note that for r = 0, 7t t> q(s) reduces to the ideal gas partial pressure 
?r(s), c.f. (3.9). 

Define the following probability measure on Q: 

We shall consider the cases < e and Ek > e individually. To do so, define a pair of 
occupation measures as follows: 

so that both L measures are positive and bounded. Now define £\ := Ai+([0,e)) and E% := 
J\A b + ([e, oo)) as two spaces of positive bounded measures on [0, e) and [e, oo) respectively. 
We equip both these spaces with the narrow topology. Then set £ := Si © £ 2 . 

Let K a 'W (for i = 1, 2) be the probability measure induced on £j by L: 
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where for brevity we will not write all the arguments, and set 



Then define the Legendre transform 

TV*(t)= SUp {tS-TT(s)} =tS*(t) -7T(s*(t)) 
s<miii{0,— t} 

where s*(t) < is the maximiser satisfying 

t = n'(s*(t)). 

Similarly set 

<q0) = sup {ts(t) - n T>Q {s(t))} = ts* TQ (t) - n T>Q {8* iQ {t)) 

s<min{0,— r} 

where s*n(t) < satisfies 

t = < Q « Q (t)). (3.19) 

As above for i — 1, 2 and every € let m^ = mi +mi be the Lebesgue decomposition 
of with respect to dz/, where mi is the singular part and mi* ^ the absolutely continuous 
part. Define [/, : Si —> [0, oo] by 

U?[m] = - [ (a- \)m^(d\), U?{m] = - [ (a- \)mf\d\). 

J[0,e) J[e,oo) 

As above, fix 

J 1 (t,r)=ir*(t)-rt + ir(r), J 2|T (t, r) = ir* )Q {t) - rt + ir T>Q (r). (3.20) 

Then by Theorem 2 of [39], the sequence of measures satisfies the Large Deviation 
Principle with constants Vi and rate function / of the form 

where 
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I°[m] = U?[m] + / Ji(p(X), a - \)dv{\), 

'[0,6) 



j« r [ m ] = [f«[ m ] + / J 2 , T (p(A), a - X)du(X). (3.21) 

./ [e,oo) 

Then we may use Varadhan's Theorem to show that the pressure of the modified Perturbed 
Mean-Field model (using the a < trick) in the thermodynamic limit can be expressed as 

p(p, r, Q, e) = p°(a, r, Q, e) + sup [G^ a [m] - J r Q [m]} 

= p°(a,T,Q,e) + sup {G^K] - IfK] + G^ Q [m 2 ] - J 2 a T [m 2 ]} 

m=mi+m2 

= p°(a, r, Q, e) + G M_Q [m*] - /?[<] - ^Kl 

for some m* e £. The Euler-Lagrange equation for the variational problem implies 

d d 
— p(p, t, Q, e) = ^rP°(« 5 T , Q, £ )+ 

^« - ^m<] - ^vKl- ^)) - ^LKl (3-22) 

of which the inner product must be zero because j- (G[m*] — J"[m*] — I^[m*]) = 0. Then 
taking the remaining term, and using (3.20) and (3.21), 



d_ 



!U m *\ = f |"-V(p(A), a - X)du(X) 

J[e,oo) ° T 



d f d 

7T* TtQ (p(X))du(X) + / 7r-7r TiQ (a - X)du(X) 



e,oo) 9t T ' J[e,oo) 

—^ Q { P {X))du{X) + ^-P?(«, Q, r, e) 

[e,oo) ^ ^ 

since n has no r dependence. Inserting this into equation (3.22) causes the free pressures to 
cancel, leaving one with 



d_ 
dr 1 



f d 

p(p,T,Q,e)= —n* Q (p(X))dv(X). 

J[e,oo) OT 
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But note that 



<«(*) = (* " <«««(*))) 57 - ^,Q«q(*)) = 



<9r 



<9r 



since t = 7r^g(s*g(t)) from (3.19), which returns the following expression for (3.22), the 
derivative of the interacting pressure with respect to r: 



where K.o satisfies t = ir' Q (s* Q (t)). 







[e,oo) 



^ 7r r,Q« Q (p(A)))rfz/(A). 



(3.23) 



From (3.18), the definition of 7r Tj g, one can can immediately calculate that 



d_ 

dr 



T=0 



-Ps* _ I ■ 



(3.24) 



Note that when r = 0, vr r= o.g reduces to the ideal gas partial pressure (3.9) for which we 
have the explicit expression: 

7r r=0i g( S ) = 7r( S ) = iln [1-e*] -1 

which is independent of Q (so we drop its index). Differentiating this with respect to s 
implies that 

t 



t 



l + t 



and thus the integrand of the non-zero term of (3.24) in terms of A is 



>PQs* (p(\)) 



-/J«S(pW) _ i \ p (X) + l 



P(A) 



p(X). 



Combining this with (3.23), 



lim yV(g,/-P f 

q=Q 



) = - 



p(/i,Q,T,s) = / I 

=0 J{s,oo) \ 



P(A) 



P(A) + 1 



p{\)dv{\) 



which goes to zero at Q — > oo by Lebesgue's Dominated Convergence Theorem. 



□ 
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Finally we note that when e ^ has a Feynman Kac representation, Condition III follows 
from (3.2) since by (3.3) 

limlim c»(q,P £ ) < limlime^ 9 [ dvAX) < C d e? M lim I \ d / 2 d\ = 

using the expression for the integrated density of states (3.5). 

However for more complex models not all these conditions are not necessarily true. In fact 
we know that the vital underlying assumption of the non-negativity of the cycle expectation 
for a positive operator (Condition I) is not satisfied in general. In the following two chapters, 
we shall consider the Infinite-Range-Hopping Bose- Hubbard model in two forms, first with a 
hard-core repulsion explicitly forbidding more than a single particle per site and second with 
a finite on-site repulsion to discourage multiple particles per site. In Appendix B, by direct 
computation we find that the positivity condition does not hold for the hard-core model, 
and therefore an alternative approach is required. 
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I.R.H Bose-Hubbard Model with 
Hard-Cores 

Summary 

In this chapter we study the relation between long cycles and Bose-Condensation in the 
Infinite-Range-Hopping Bose-Hubbard Model with a hard core interaction [15]. Unfortunately 
the technique derived in Chapter 3 cannot be applied to this model ( as Condition I fails to 
hold, see Appendix B). Nonetheless we succeed in calculating the density of particles on 
long cycles in the canonical thermodynamic limit by the following argument. We express the 
density of particles on cycles of length q forn particles in terms of q distinguishable particles 
coupled to a field of n — q bosons. We can prove that in the thermodynamic limit we can 
neglect the hopping of the q particles so that bosons have to avoid each other and the fixed 
positions of the distinguishable particles. This is equivalent to a reduction of the lattice by q 
sites. Moreover the q particles are on a cycle of length q. It then is shown that in fact only 
the single-cycle density contributes which means that in the thermodynamic limit the sum 
of the long cycle densities is the particle density less the one-cycle contribution, which can 
be calculated. We find that the existence of a non-zero long cycle density coincides with the 
occurrence of Bose-Einstein condensation but this density is not equal to that of the Bose 
condensate. 
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4.1 The Model and Previously Derived Results 

As shown in Chapter 3, the density of long cycles is equal to the density of the Bose-Einstein 
condensate in the cases of the Bose gas, the mean-field model and the perturbed mean-filed 
model. We wish to check if this is the case for the selected model: the Infinite- Range- Hopping 
Bose-Hubbard Model with hard-cores, also simply referred to as the Hard-Core Boson model. 

The Bose-Hubbard Hamiltonian is given by 

h^ = j « - «pk - <> + A E n ^ n * - x ) ( 41 ) 

\x-y\=l 

where A v is a lattice of V sites, a* and a x are the Bose creation and annihilation operators 
satisfying the usual commutation relations [a x , a*] = 5 X)V and n x = a* x a x . The first term with 
J > is the kinetic energy operator and the second term with A > describes a repulsive 
interaction, as it discourages the presence of more than one particle at each site. This model 
was originally introduced by Fisher et al. [43]. 

The infinite-range-hopping version of this model has Hamiltonian 

RlR = ^7 E (°i - a i^ a * - s) + A E n ^ - x )- ( 4 - 2 ) 

This is in fact a mean-field version of (4.1) but in terms of the kinetic energy rather than 
the interaction. In particular as with all mean-field models, the lattice structure is irrelevant 
(due to the infinite-range hopping) and there is no dependence on dimensionality, so we can 
take A v = {1, 2, 3, ... , V}. The non-zero temperature properties of this model have been 
studied by Bru and Dorlas [17] and by Adams and Dorlas [1]. Dorlas, Pastur and Zagrebnov 
[30] considered the model in the presence of an additional random potential. However this 
chapter will consider a special case of (4.2), introduced by Toth [80] where A = +oo, that 
is with complete single-site exclusion (hard-core). The properties of this hard-core model 
in the canonical ensemble were first obtained by Toth using probabilistic methods. Later 
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Penrose [66] and Kirson [57] obtained equivalent results. In the grand-canonical ensemble 
the model is equivalent to the strong-coupling BCS model (see for example Angelescu [5]). 

We recall the thermodynamic properties of the model in the canonical ensemble as given by 
Penrose (see Appendix A for a summary). 



For p E (0,1), let 



ln(— P -\ 



9(P)={ 1 - 2P V 9 ' (4-3) 
2 if p = 1/2. 

For each (3 > 2 the equation (3 = g(p) has a unique solution in (0, 1/2] denoted by pf r (see 
Figure 4.1). We define p(? r := 1/2 for (3 < 2. The inverse temperature (3 will be fixed for the 
remainder of this chapter, and so in maintaining the notational custom of this document, 
we shall only write the (3 subscripts when necessary. 

Theorem 4.1 (Penrose [10], Theorem 1) The free energy per site at inverse temperature [3 
as a function of the particle density p e [0, 1] , f{p), is given by 



f(p) 



p+-(plnp+(l-p)ln(l-p)) if pe [0,p cr ]U[l-p cr ,l] 

p 2 + p cr (l - p cr ) + - (p cr lnp cr + (1 - p cr )ln(l - p CT )) if p G [p cr , 1 - p cr ]. 



The density of particles in the ground state in the thermodynamic limit is given by 



n/v=p x,y&A v 

where (■) denotes the canonical expectation for n particles, c.f. (2.10). Penrose showed 
that for certain values of p and j3, Bose-Einstein condensation occurs, that is ^ond > ^- The 
Bose-condensate density is given in the following theorem. 

Theorem 4.2 (Penrose [10], Theorem 2) 

The Bose-condensate density, Q p cond at inverse temperature (3 as a function of the particle 
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g(p)k 




density p G [0,1], is given by 



Figure 4.1: Definition of p c 



Second 



if p e [0,p cr ] u [i-p cr ,i], 

(P - Pcr)(l - P ~ Per) if P 6 [Per, 1 - Per]- 



(4.4) 



We note that both /(p) — p and the condensate density Q p conA are symmetric about p = 1/2. 
This can easily seen by interchanging particles and holes. The boson states being symmetric 
can be labelled unambiguously by the sites they occupy but equivalently they can be labelled 
by the sites they do not occupy (holes). 

It is also worth noting that this model is equivalent to an XY model on a complete graph. 
The mean field theory of the XY model was calculated by Matsubara and Matsuda [65] 
who found an formula for the critical temperature equivalent to p cr , and Zitsel [101] gave 
a free-energy formula equivalent to /(p), and a density-temperature graph essentially the 
same as Figure 4.1. 

We wish to study the cycle statistics of the Infinite-Range-Hopping Bose-Hubbard model 
with hard cores and compare Q p conA with the density of long cycles. 
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Some Notation 



Before proceeding we need to define the n-particle Hamiltonian more carefully. The single 
particle Hamiltonian has eigenvalues and 1. Denote the basis of the single particle Hilbert 
space of V sites, 1-L V := C v , by {ej}J^ 1 where e; = (0, . . . , 0, 1, 0, . . . , 0) with 1 in the 2-th 
place. A Hamiltonian for this model is the following: 

h v \<j>) = \<j>) - (g\<f>)\g) 

or simply expressed as 

h v = I- P v 

where P v = (g v \ ■ )g v is the projection onto the ground state g v defined by the eigenvector 

with V — 1 other possible (excited) states Cj such that Cj _L g ( c i|g) = 0. Then P v is 
given by 



Thus h v is the orthogonal projection onto the subspace orthogonal to g v . With this notation 
we can define the non-interacting n-particle Hamiltonian hy acting on the unsymmetrised 
Hilbert space as: 

h (n) = jin) _ pin) 

= n-P v ®I®---(g)I-I®P v (g)I(g)---®I I <g> I <g) • • • <g> P v . 

For bosons we have to consider the symmetric subspace of %[™\ which we have denoted 
When hy is restricted to 7iy + , we obtain the usual free Hamiltonian kinetic energy: 

^yYl ( a l- a *y)( a x- a y)- 

x,y£Ay 
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We introduce the hard-core interaction by constraining the Hamiltonian to a subspace of T-Ly 
where states with more than one particle per site are discarded. This is done by applying 
the following hard core projection to 

{0 if e ik = Bi , for some k ^ k', 

e h <g e i2 <g ■ • • <g e in otherwise. 

We shall call H^y '■= V^Ti^ the unsymmetrised hard-core n-particle space and 
'H h nv+ := T*n C %v + t ne symmetrised hard-core n-particle space. Note that as [U w , V* c ] = 
for all 7r G S n , V^ commutes with the symmetrisation and so 1-L n c v + = 0+^ 7~L n c v - 

The hard-core n-particle Hamiltonian is then the following operator 

H% ■= K C h { v ] K c (4.5) 

acting on the hard-core n-particle space H^ v - Therefore the Hamiltonian for the infinite- 
range Bose-Hubbard model with hard-core is (4.5) acting on the symmetric hard-core n- 
particle space %^ c v+ . 



4.2 Applying Cycle Statistics to this Model 

We wish to apply the framework established in Chapter 2 to this model. The canonical 
partition function for this model is 



Z v (n) = trace W hc 



n,V,+ 



e -P H n,V 



1 



— trace -i/hc 

fp Z— / n n,\ 



7r65„ 



Recall that the validity of this framework rests on (2.2), i.e. the condition that each term in 
the above summation is non-negative: 



trace u^c 

H n,V 



> 
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for all 7r £ S n . But from the random walk formulation (see for example [80]), one can see 
that the kernel of e~^ H ^ v is positive and therefore the left-hand side is always positive. 

Following the procedure of the proof of Theorem 2.1 with the following replacements: 

hi -> h% v n { : ] v^u { : ] := U% 

Z A (n) Z v (n) ® V* c {%f ® U^ q) ) 

one may derive the the following more expression for the g-cycle density: 



C " (9) = W traC6 ^ 



where H^ v := V^(U^® H%?>) 
Then as before we have that the density of particles on cycles of infinite length is defined by 

oo 

<n g = lim Km Yl c v( q } 

n/v=p q=Q+l 

(c.f. Definition 1). We wish to compare this quantity with the condensate density Q p conA - 
By defining the density of short cycles as 

Q 

0short = P - <ong = n lim 1™ XI C "^' 

n/v=p 9=1 



then observe that we can interchange the q summation and the thermodynamic limit to 
obtain the expression 

Q 

£short = n lim Yl I™ 
9=1 n/v=p 

which will make it much easier to calculate later. 
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4.3 Neglecting the Hop of the g-unsymmetrised Particles 



By using cycle statistics we can split our symmetric hard-core Hilbert space H^v + 
a tensor product of two spaces, an unsymmetrised g-particle space T-Ly and a symmetric 
n — q particle space Hy \ , with the hard-core projection applied. However in doing so 
we explicitly isolate a finite quantity of distinguishable particles. In the thermodynamic 
limit, the total number of particles goes to infinity, so one queries whether the hopping 
contributions of the q distinguishable particles are significant in the limit. 

To isolate the hopping contribution of these q unsymmetrised particles, we need to establish 
some notation to clearly segment the action of an operator onto the two Hilbert spaces Ti^ 
and H { yl q) . Write 



A (q) - = A (q) g, jin-q) &nd A (n- q ) ._ j(g) ^ A (n-g) 



/he 



for any operator A on T-L v . With this form we can express the Hamiltonian (4.5) on T-L qn 
by the following: 

Then the operator P v q ^ = V^P V represents the hopping contribution of the q distin- 
guishable particles to the Hamiltonian (4.5). 

Now let us define the following modified Hamiltonian 

H&v = Pn (" " Pv l ~ Q) ) Pn (4-6) 

so that 

trhc _ frhc pil) 

n,v n q,n,v V ' 

This Hamiltonian considers n particles on a lattice of V sites, with q particles in fixed 
locations, and the remaining n — q particles hopping from site to site with the hard-core 
constraint (a particle must avoid landing on any other particle). One may observe that this 
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ensemble is similar to a system of n — q hopping particles on V — q sites with hard-cores, a 
fact that will be made rigourous shortly. 

Now we wish to estimate the effect of neglecting the action of the P v term in the unsym- 
metrised space. Define the modified cycle density by 

= ^^trace^c 

We want to show that the g-unsymmetrised particles' hopping contribution may be neglected 
in the limit. This is expressed in the following theorem 



Theorem 4.3 

lim Cy{q) = lim c"(g). 

n,V— too n,V — >oo 

n/v=p n/v=p 

Proof of this result requires two lemmas. In Lemma 4.1 we find an upper bound for the 
difference of the two cycle densities for fixed q, which involves a ratio of partition functions, 
one of n particles on V sites, the other oin — q particles on V — q sites with a scaled hopping 
parameter. Secondly in Lemma 4.2 we determine the limiting behaviour of this ratio. 

A little more notation before proceeding, define a scaled partition function by 

Z v (\,n) = trace ^hc^ 

where the Hamiltonian introduces a weighting parameter A to the hopping action: 

Aw = (n - \Pi n) ) V^. (4.7) 



The first step is determining the following estimate. 



Lemma 4.1 

\c«(q)-c«(q)\< {1 -^ 9)Zv -^ n - q) 
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To prove this Lemma we intend to obtain an upper bound for 



U q e 



U q e 



^ q,n,V 



However before proceeding with a proof we first shall introduce some notation and make 
some remarks. 



Some Preliminaries 



Let 9 ' be the family of sets of n — q distinct points of A v . For k = {k\, k 2 , . . . , k n _ q } £ 

(n-q) 

|k) := <7+~ ff (e fcl ® efca ® • ■ ■ ® e kn _ 



A^+ q) let 



Then {|k) | k £ Ay\ ^} is an orthonormal basis for H^_ q v + ■= V^ c _ q l-iy l + . 

Similarly let Ay be the set of ordered g-tuples of distinct indices of A v and for i = 
(21,22, ...,i q )e A iq) let 

|i) := e h ® e i2 <g> • • • <g) e v 
Then {|i) | i £ A^} is an orthonormal basis for "Hjl c v := V q c Hy\ 

If k £ A[ n + ^ and i £ A^ we shall write k ~ i if {fci, fc 2 , . . . , k n - q } Pi {ii, 12, ■ ■ ■ , i q } = and 
we shall use the notation 

|i;k) := |i) <g> |k). 

Then a basis for V^ qnv may be formed by taking the tensor product of the bases of %^L q v + 
and Hq V where we disallow particles from appearing in both spaces simultaneously. Thus 

he 

q,n,V 



the set {|i; k) I k £ A^ q \ i e Av\ k ~ i} is an orthonormal basis for % 



We shall need also the following facts. For simplicity we shall write H and P for H q c n v (c.f. 
equation (4.6)) and Py^ respectively. 
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Let , p} n ~ q> be the projection of W£_ q <, v + onto a space with none of the n — q particles at 
the points ii, 12, ■ ■ ■ , i q (so there are V — q available sites for n — q particles) and not more 
than one particle at any site. Then we state the following useful remarks: 

Remark 1 For i ~ k, if s > 0, then 

e-^ s |i; k) = |i; e-^ is k)e-^ 9S (4.8) 
where H 1 = v[ n ~ q \{n - q) - p^)p^ . 

Proof: This can be seen as follows: For i ~ k, 

H\i;k)=V h n c (n-Pir q) )V h n c \i;k) 

= qV h n c \i; k) + 7>£=|i; ((„ -q)- P^)k> 
= g|i; k) + |i; v[ n ~ q) {{n - q) - P^ q) )k) 
= g|i;k) + |i;iTk). 

□ 

Remark 2 For i ~ k, 

n-q 

H l \k) = (n - q) |k> - - Y, E I *&> ■ ■ ■ > ^ ^' ■ ■ ■ » 

Z=l j£iUk\{fc ; } 

where the hat symbol implies the term is removed from the sequence. 



Proof: This statement follows from the definition of H 1 and the properties of the ket. Note 
that the I summation selects each of the n — q particles in turn, and the j summation is 
over all the possible landing-sites for each particle (i.e. that do not violate the hard-core 
condition). □ 
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Note that from (4.7), for k G T-Ly_g + we have 

n-q V-q 

Hx%-q^q\V = (n- q )\k)--—T E \(h,k 2 ,...,k l ,j,...,K- g )). 

y 1 = 1 3=1 

j^k\{fc ; } 

Thus if 1 is unitarily equivalent to y v n _ qv _ q and 

= ^v- g (np, n-g) . 



trace ■l/hc 



'p( n -l) e -PH l <p(n-q) 



Remark 3 For s,«GR 



^(n-gJg-sH'pCn-g)^ 01 _ <p{n-q) e ~saH l <p{n-q) 

Proof: Let A = Vie~ sH 'Vi, B = e~ sn \ We shall show that for any / continuous, f(A) = 
Vif{B)Vi- These are self-adjoint operators, so using the Spectral theorem, there exists 
unique spectral measures fj,^(dX) and (jJp JdX) respectively so that 

(<f>\A a <f>) = (VMB^) => /xJ(dA) = ^(d\). 



Therefore we may write 



WW® = I /(A)^(dA) = J f(X)$ BVi (d\) = J f{X)4 i<f> {dX)^{V i( j>\f{B)V4)- 



By setting f(x) = x a we retrieve our result. □ 



Proof of Lemma 4.1 



To proceed we expand 



trace o/hc 

1,n,V 
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-l/hc 

q,n,V 



U q e 
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in a Dyson series (see [73; 74]) in powers of P. If m > 1, the m th term is 



X m :=r [d Sl [ds 2 
Jo Jo 



(I/Sim Xxbjzq -^hi 



q,n, V 



e -PH(l- Sl ) p e -fSH( Sl -S2)p , 



pQ—/3H(s m -i—s m )pg—f3Hs m ]J 



Recall that P := V^P^V^ where 



P® = P V ®I 



I + I®P V ®I 



I+--- + I®---®I®P V 



has g terms, so in the above trace we have m instances of this form. Let = I 



P 



r place 



® /, and let P r = V^P^V^ 



Then we can write 



where 



x m = r J2 E ■ ■ • E x ™( r i' r 2, . . . ,r 9 ) 



(4.9) 



ri=lr 2 =l r m =l 



X m (ri,r 2 , • 

x trace 



1 />Sl 

I ds\ ds 2 . . . I ds 
o Jo Jo 



m — 1 

m 



he 



-PH(1-s{) p -/3H(s 1 -s 2 ) p ...p „-/3H(s m - 1 -Sm) p -&Hs m TT 



In terms of the basis of "H^ c n v we may write the expression for X m (ri, r 2 , ■ ■ ■ , r q ) as 



X m (ri, r 2 , . . . , r 3 ) = /dsi ds 2 . . . / ds m ^ E " ' E 



E 

k ,... ,k m i°~k° i m ~k" 



>1 /-Si 

dsi / (is 
'o Jo Jo 

(i°; kV^ 1 "' ^ 1 ; k 1 )^ 1 ; kV^^^li 2 ; k 2 ) ■ • • 

. . . ^-ljg-^C^-i-^p^ijm. k ™)(i™ ; k m | e -^ Sm t/ g |i°; k°) (4.10) 



where it is understood that the i summations are over Ay and the k summations are over 



A 



(n-q) 
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Note that for i ~ k 

P r |i; k) = 1 l(*i> • • ■ ^r, I,---, i g )\ k) (4.11) 

1=1. ..v^ 
10C', l^i\...i T ...i q 

where again the hat symbol implies that the term is removed from the sequence. 

Consider one of the inner products in the expression (4.10) for X m (r 1; r 2 , . . . ,r q ), using (4.8) 
and (4.11) above. For i ~ k and j ~ k': 

(i;k|e^P r |j;k') = ^ £ (i | ( Jl; • • • , Jr , • • • , j 9 ))(k | e^|k') 

1=1.. .v 
l(£k';l^j 1 ,...j~ r ,...j q 

= — E S iljl ...6 irjr 6 irl ...6 iqjq {k\e-^ H \k'). 

1=1. ..v^ 

l(£\i';l^j 1 ,...,jr,--;jq 

In summing over I we replace / by i r and the result is non-zero only if i r ^ k' and i r ^ 
ji, . . . , j r , . . . ,j q . However this last condition is not necessary because if i r = j s (s ^ r) then 
j s 7^ i s and we get zero. Also if for some s / r, i s e k' then once again j s ^ i s . We can 
therefore replace the condition i r ^ k' by i ~ k'. Using I for the indicator function, we have 

~ _ e -Pqs . . 

(i;k|e-^P r |j;k') = . . . ^ <k|e-^" |k'>^(i~k') 

= . . . C • • • S iqjq (k\vt q) e-" sHi vt q) \V)- 

Now if we sum over j ~ k', with i ~ k and for a fixed r: 

^(i;k|e-^P r |j;k')(j;k'| = ^<k|^-V^^- 9) |k'} 

X ((ii, . . . , % r -\i jr, Vflj • • ■ j iq)'i k |- 

Jr=l—V 
j r <£U.'Ui\{i r } 

It is convenient to define the operation [r, x](i) which inserts the value of x in the r th 
position of i instead of i r . So for example taking the ordered triplet i = (5,4,1), then 
[2, 8] (i) = (5,8,1). For brevity we shall denote the composition of these operators as 
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[r k ,x k ; . . . ■,r 2 ,x 2 ]r 1 ,x 1 ] := [r k ,x k ] o ■•• o [r 2 ,x 2 ] o [ri,Xi]. 

At a very low level, this operation represents the individual hopping of the particles. To 
prevent violation of the hard-core condition, one needs to take care when constructing the 
set of all possible landing sites. 

In this notation the final term in the above expression may be rewritten as ([r, j r ](i); k'|. 
Performing two summations for fixed r\ and r 2 we get: 

J2 (i^^le-^Pnli^k^^jk^e-^P^Ii^k 2 )^ 2 ;^! 
ii^k 1 i 2 ~k 2 

e -Pq(s+t) 



E (kVioe-^ T^k 1 ) 



V 2 



x (k l \V [ri , lhm e -^ l4l](i0) ^^jciojlk 2 ) <[r l5 4 X ; r 2 , ^ 2 ](i°); k 2 | 



Y, E < k °i <^"" i kl > 

x (k 1 ! P [rii4i](i0) e -^ [ri4lKi0) P [n , 4i](i0) |k 2 )<[r 2 ,^;r ll tj l ](i );k a | 

due to the fact that P^jij^o^k 1 ) = if i\ x G k 1 . We may apply this to all inner product 
terms of (4.10) except the final one. Note we sum over the V sites of the lattice, with certain 
points excluded in each case. 

For the final inner product of (4.10) we obtain: 
([r ro ,i^; ... ;r 2 ,i^;ri,iy(i );k TO |e-^™%,|i°;k > 



• <[r m , C; • • • ; r 2l £;»•!, ^J(i°); k-| e -^ [rm '^- ^'^Vji ; k o } 



e -^ m(k m |e -^'-— S-^) |k0) ([rm , C; ... ;r2 ,^. ri?z i J(i o )| ^ i o ) 

o-^m/L-ml-T, — /JSm-fi'' 1 "" 1 '*"' 71 ' "' ;r 2' i r 2 ' r l' i r 1 ]( i0 ) 

e ( k lHr m ,i™ m ;...;r 2 ,i2 2 ;r 1 ,il 1 ](iO)e HwP m ! ••• in.f^KiO) I k ) 



x ([r m ,C; ... ;r 2) tJ,;r 1 ,iJ 1 ](i )|U g i > 
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Applying this to the whole tracial expression of (4.10) we obtain 



X m (r 1 ,r 2 , ■■■,r q ) 



C 

y m 



E E E E E 

(k^Pioe^^-^^Violk 1 ) 

-0( B i- S2 )tf [ri,i M |O >, 



^I^In.i^KiOje" 



^[n,e ](i°)l k ) 



i 1(i o )e -^---)^ lr2 ' V2;ri ' VlKl V 



r2,i?„;ri,i; 



,](i»)|k 3 ) 



(k m |P [r 
([r m ,i 



/9 am H lrmi *^ ; ' ■' ;r2 ' l? 2 ;n '^l 1 ( ' 0) <r> |T_0\ 

Hr m ,»™ ;...;r 2 ,ip, jri,** ](i°)l k ) 



i? 2 ;ri,ii 1 ](i°)e 



,...;r 2l ^;r 1 ,ij 1 ](i )|^i > 



-#7 



V 



EE E 



E 



{[r m ,C; ...;r 2 ,^;r 1 ,z r 1 1 ](i°)|t/ g i ) 



n (n-q),V,+ 



r-2 1 ' J-' "nJ /rf 

.-.;() 



p i oe-^ 1 - sl ^ 1 > i o^i i(i o ) e-^ sl - S2 ^ Vl(i) 

2 .1 1,-iO 



^1 (iO) 



Pi 



From the Holder inequality (see Manjegani [63]), for finite dimensional non-negative matrices 
Ax, A 2 , . . . , A m+ i we have the inequality 



771+1 



[trace (AiA 2 . . . A m+1 ) | < trace j^x^ . . . A m+1 | < 1 [ (trace A^ 



j_ 



k=i 



where YJk=i h = x > Pi > 



Pk 



Set pi 



l-si 



P2 



S1-S2 



j • • • ) Pm 



p m+ i = — . Taking the modulus of the above 
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trace 



n (n-q),V,+ 



-(3H i0 (l~ Sl ) 



Vioe 

'^[fm,^ m i».;r2i'?/l,<y( i() ) e 



VioV [riAim e pH 1 ri (S1 S2) V [riA i im ■ ■ ■ 



I fr™.j?P :...:ro.i? :ri .ii lCi°1" ' [r m ,i™ m \ ... ;r , 2,*r 2 ' r l'*r 1 ](i°) 



< trace 



"(71-9),^,+ 



Pioe"^ 1 7>, 



trace 



trace 



5[n,irJ(i U ) 



V 



81- 



^Vm,iS_ ; ••• ;^2,iL;n,iri l( i0 ) e 



;r-2,C;ri,ji ](i°) 



But since the trace is independent of the 1/ — g sites ](i°), [r 2 , ri, ^J(i°), . . . 

. . . , [r m , ; . . . ; r 2 , 2^ 2 ; r 1; ^J(i )}, therefore using Remark 3, the product of all the trace 
terms above is equal to 



tjI*clC6 



-(n-q),V,+ 

with 1 = {V — q + 1, V — q + 2, . . . , V} and from Remark 2, 



tjI*£lC6 



( n -q),V,+ 



ne-^Pi] =Z v _ q (^,n-q) 



The integrand is now independent of s±, . . . ,s m , so we may integrate to obtain a factor of 
\jm\. What remain is to consider the unpleasant looking sum 



EE E 



E 



([r m ,iZ; ... ;r a) ^,;r 1 ,if l ](i )|^i ). (4.12) 



t [rm-i ,C m -_\ i ... j^,^ jn,^ ] } 



To analyse this, first fix the numbers r 1; r 2 , . . . , r m . If {r 1; r 2 , . . . , r m } ^ {1,2, ... ,q}, then 
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the ket |[r m ,i™; ;r 2 ,^ 2 ;ri,^ 1 ](i )} is of the form 

jli • • • , jni , • • • ' ^n 2 ' jn 2 +li • • • i Jnzi ^n 3 +l5 • • • ) ^7141 jVi4+l> } 

where {ni,n 2 , . . . } is a non-empty ordered set of distinct integers between and q. This 
state is clearly orthogonal to U q i° for any q. Note that this situation does not arise if q = 1. 
Note also that this is always the case if m < q. 

We may bound the remaining sum corresponding to terms for which {7*1, r 2 , . . . , r m } = 
{1,2,...,?} by 

^E E E ••• E (K,c;---;^v 2 2 ;n,^](i o )i^i°). 

:0 ,-1 —1 ,;2 —i im —\ 

v ' 

where [r m , ; ■ ■ ■ ; ^2, ; fl , ] 
has distinct indices 

Observe that in this case | [r m , i™ m \ ... ; r 2 , if ; 7*1, ^J(i )) is independent of i° so we may take 
it to be 

••• ; r 2,^ 2 ;n,^J(s )) 

where s° = (1, 2, 3, . . . , q). Then we can interchange the i° summation with the others, and 
for each choice of , if 2 , . . . , ijf 1 there exists only one possible i° G such that 

([r m ,C m ;...;r 2 ,4;r 1 ,4 1 ](s°)|C//)^0. 

So we may conclude that 

EE E 

(K,C ; ••• ^il^rUj?)\u q ?) <v m . 

^^[r™-i,C m -_ 1 1 ;...;r 2 ,i2 2 ;r 1 ,ii 1 ](i«)\{C m - 1 } 

Applying this, we see that the modulus of the integrated m th term of the Dyson series may 



76 



CHAPTER 4. I.R.H BOSE-HUBBARD MODEL WITH HARD-CORES 



bounded above by 



\X m (ri,r 2 , ...,r r 



< 



V m m\ 



*-.(V.»-«)E E E • • 

^ m £ [r m - 1 ,i^T_\ ; • • • ;r 2 ,ip 2 ;r i ,i\ ] (i°)\{CT 1 } 



m! 



-3) 



Hence the modulus of (4.9), the m term of the Dyson series may be bounded above by 



\X m \ < P m J2 ■ ■ ■ E \ X ™(ri, r 2: . . . , r m ) 

ri=l r„ = l 



< e^Z v _ q (E=2, n -g) 



m! 



Noting that the zeroth term of the Dyson series is 



X n = trace 



nj he 



we may re-sum the series to obtain 



trace 



H 



he 

q,n,V 



U q e 



<e-^Z v . q {^,n-q) £ 



m=l 



ml 



Thus 



c v {q)-c v {q)\ 



V 



U n e-?<v 



Z v (n) 



< g z ^-g n-q) y-v g p 



m=l 



m! 



,/3g 



Zy_ g n-g) 



Z v (n) 



(4.13) 



Obviously we want this to tend to zero in the limit, but this is far from clear. The next 
section will consider the behaviour of the troublesome ratio of partition functions in the 
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thermodynamic limit. 



4.4 The Variational Argument 



For this section we shall reinsert the [3 dependence in all quantities where required. 

In his paper on the Hard-Core boson model [66], Penrose gave an explicit expression for 

Z v (n,f3): 

min(n,V—n) 

Z v {n,p)= <r,n,V,P), 

r=0 

where 

z(r,n,V,f3)= ( j^^f j (^) ex p{"f [Vr - r 2 + r + n 2 - n]^. 

He also proved that if h v : [0, min(p, 1 — p)] — > R converges uniformly in [0, min(p, 1 — p)] as 
V — > oo to a continuous function h : [0, min(p, 1 — p)] — > R, then 



_^ mm(n,V—n) 

n ,v^ Z v (n,B) ^— ' v 



%), if pe [o,pf r ], 

% cr ), if pG^lVi ( 414 ) 
/i(l-p), if p<=[l-p£,l]. 



We shall use this expansion and his Large-Deviation-like result to consider the limiting 
behaviour of the ratio of partition functions in line (4.13), to obtain the following: 

Lemma 4.2 



, Z v _ q (^,n-q,(3) \p q ^ q if pe [0,/4]U [1-/4,1], 

hm — 



n/v= P *W> p > I eftrd+z-A) if p e [^ 1 



Z v (n,/3) - p 13 ). 
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Proof: Recall that we have 



Zv- a (n — q, 8) = trace -uhc 

n-q,V-q,+ 



while 



-BH hc 

' n — q,V—q 



6 ^ ^ trclCG ay he 

rL n-q 1 V-q,+ 



/□-phc p n — <?«phc 
i^' n-q r V-q' n-q 



(4.15) 



(-^f 5 , n- q,8) = trace 



-l/hc 

rl -n.-q,V — q,+ 



-/3(n-g) 



g/^V y n-g^V-g' n-q 



(4.16) 



Comparison of (4.15) and (4.16) yields an equivalence between the partition functions 
Z v _ q (V, n-q, 8)= e-fV-^) Z v _ q (n - q, 8{^)) 

and thus we have to analyse the following ratio: 

e-W"-g>Z v _ g (n- g>/ 9(V)) 
Z v (n,/3) 

We wish to express this in the form of the lefthand side of (4.14). Using Penrose's explicit 
form for the partition function, we have the following: 



Zvr_ ff (n-g,0(V))= E 



= ( V-g-2r + l } (V - q 

r 



r=0 



V — q — r + 1 



exp | — ^ [r(V — q) — r 2 + r + (n — q) 2 — (n — g)] | 



For the case p > |, for large K, n-g>F-nwe must sum from zero to V — n and a 
straightforward calculation then gives 



V-n 



e H v 



^77( n -9) 



r=0 
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where 



h v (x) 



l-2x-(q-l)/V \ ( 1-x + l/V \ 
l-2x + l/V J \l-x-(q-l)/V J 



3-1 

X 

s=0 



Therefore 

= lim /i v (x) = (1 - x) 9 e qf}{x+p) . 

V— >oo 

It is clear that the convergence is uniform since h v (x) is a product of terms each of which 
converges uniformly on [0, 1 — p] for p > |. Thus 



lim 



e-m^Zv^ (n - q, B{^)) (1 - PcV) 9 ^ {p " +p) if P G (1/2, 1 - 



n/v=p V{ ' P) Xffl^l if p 6 (1- P f r ,l]. 



Note that using the relation 

1 - p* 9 

rcr 



/3 = ^r^ln 



1 2p C r \ Per 

we get 

(1 - p^) 9 e q ^ + ri = {(f ct )i 
and therefore we have proved Lemma 4.2 for p > |. 

For the case p < | we have that n — g < V — n, the sum for e~^v^ n ~ q ^Zy^ q (n — q, 3(^rr-)) 
is up to n — q, and therefore we need to shift the index by q to get it into the required form. 
After shifting we get 



e -0f(n- g ) 



X V(V-l)(V-2)...(V-q+l) eXP \V [V+n rJ 
Note that summand is zero if we put r = 0, . . . , q — 1. Thus we may sum from zero to n to 
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get as before 

n 

e-W^Zv-, (n - q, P(^)) = M£) z(r, n, V, 0) 

r=0 

where this time 

fl-2x + (q + l)/V\ q T i fx-s/V\ fo M n . < 



so that 



h(x) — lim h v (x) = x q exp{/3g(l + p — x)}. 

71, V — >CtO 

n/v=p 



Convergence is again uniform on [0, p] for p < | and therefore 

e -^(- 3 )Z v _ g (n- g ,/3(^)) fp 9 ^ if pG[0yj, 

lim — - = < 

n/v= P Zv ^ /3) [ (p^y e ^ P -pi) if p e K , 1/2), 

proving Lemma 4.2 for p < |. The case p = | is more delicate because the first term in 
(4.17) does not converge uniformly. We can write (taking V = 2n) 

h 2n (r/2n) = h 2n (r/2n) + — - h 2n (r/2n) 

where 

3-1 



Clearly h,2 n {x) converges uniformly on [0, 1/2] and therefore 



lim 7 T ffl E^) «(r,n,2n,)9) = (p^e^ 3 / 2 "A 



r=0 



We thus have to show that 



11111 7 ) R\ E o/ 2 " (2 H 1 *( r > n ' 2n ' ^ = °" 

n->oo Z 2n (n, p) 2(n - r) + 1 
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Since h 2 n{x) is bounded, by C say, 

lim ^— — V fea "(2n) ( 2n, /3) < lim = 0. 

«->oo Z 2n (n,P) ^ 2(ra-r) + l V ' ' '^"^oo 2nV4 



On the other hand one can prove that for n — 2n 1 ^ 2 < r < n — n 1//2 and r' > n — n 1//4 

In z(r, n, 2n, (3) — In z(r', n, 2n, /3) > — Inn 

8 

for n large, so that z(r' ,n,2n, /3) / z(r,n,2n, /3) < 1. Therefore 

lim - i - V , M t } z(r,n,2n,/3)< lim O r - w "^ 

n^ooZ 2n (n,/3) r> ^ i/4 2(n-r) + l »™ ^ z(r,n,2n,/3) 

n— 2n 1 / 2 <r<n— n 1 / 2 

^ max^(r, n, 2n, /3) 

< lim ——- — r ~ n — < lim ——- = 0. 

moo n 1 ' 4 min z(r,n,2n,p) n-^oo rr/ 4 

n— 2n 1 / 2 <r<n— n 1 / 2 



□ 



4.5 The Single-Cycle Contribution 

So far we have shown that in the thermodynamic limit we can neglect the hopping of the q 
particles so that bosons have to avoid each other and the fixed positions of the distinguishable 
particles. This is equivalent to a reduction of the lattice by q sites. Moreover the q particles 
are on a cycle of length q. 

Now let us consider these q particles in terms of Brownian motion. These q particles are 
distinguishable and are all members of the same g-cycle. For q > 1, this means for example, 
that the position of the second particle at the beginning of its path is same as the position 
of the first particle at the end of its path. But since they do not hop this is impossible by 
the hard core condition. Therefore among the short cycles only the cycle of unit length can 
contribute. 
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Since the sum of all the cycle densities gives the particle density, this means that in the 
thermodynamic limit the sum of the long cycle densities is the particle density less the 
one-cycle contribution, which we can calculate. 

These ideas are mode more concrete in the following lemma. 



Lemma 4.3 c"(g) = if q > 1 and 



lim c?(i; 

n,V — >oo 

n/v=p 



p if p e [o,p cr ] u [l - p cr , l], 

p ci e^~^ if P e [p CT ,l-Pcr]. 



Proof: Recall that 



Cvio) — ^ , x 77 trace v he 
vKHJ Z v (n)V n ^ v 



Considering the trace over % qn v , expanding it in terms of its basis {|i; k)} and using Remark 



1 above, where i ~ k 



trace ■x/hc 

n q,n,V 



k i~k 



E £ W; ki 6 "^ 1 - k > = e '" q E E wi^i^w- 



k i~k 



k i~k 



For 5 > 1, an element of the basis of the unsymmetrised g-space may be written as an 
ordered g-tuple i = (ii, i 2 , ■ ■ ■ , i q ) where the V s are ah distinct. Then we may write 



(U q i\i) = (U q (e h g> e i2 (8) • • • <g> e i? ) | <g> e i2 <g> • • • ® e^) 

■ ■ ■ (g) e { ® e h \ e h (g> e i2 (g> ■ ■ ■ (g) a ) = 0. 



Hence Cy(q) is non-zero only if q = 1. 
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For the second statement, note that we may re-express c"(l) as follows: 



;trace , 
v 



Z v (n)V~^ vh ^ n v® H v: + 1) ) 



£]>>|e-nk> 



(n)V 



Zv(n)V^ 



E tracer 



.0 Zy_i(^,n- 1) 



and the result follows from Lemma 4.2. 



□ 



4.6 Conclusion 

From Lemmas 4.1 and 4.2 we found that 

lim Cy(q) = lim c™(q). 

n,V— >oo n,V— >oo 

n/v=p n/v=p 

which implies that in we may neglect the hopping of the q distinguishable particles in the 
thermodynamic limit. Therefore the density of short cycles for the Infinite- Range- Hopping 
Bose-Hubbard Model may be expressed in terms of the modified cycle density, where those 
q particles are fixed: 



£short = n lim Yl I™ c v^) 



2 — >00 n,V— »oo 

q=Q+l n/v=p 

Since by Lemma 4.3 only cycles of unit length contribute to the summed cycle density, it 
follows that 

\p if P e [o, Pcr ]u[i-p cr ,i], 

£& ort = jjm c£(l) = I 

n/v=p [pcr^ {P " Cr) If P € [p cr , 1 - p cr ] 

and therefore allows us to state the main result of this chapter: 
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Theorem 4.4 The expected density of particles on cycles of infinite length, £>f ong , at inverse 
temperature (3 as a function of the particle density p e [0, 1], is given by 



In comparing this result with the expression for the density of the condensate (4.4), we may 
form the following conclusions (see Figure 4.6): 



• £kmg * s n °t symmetric with respect to p = 1/2. As mentioned above the symmetry 
of the model about p = 1/2 is due to the particle-hole symmetry. But the equivalent 
labelling of states by sets of occupied or unoccupied sites (particles and holes) 
cannot be used for distinguishable particles. As we derived in Chapter 2, the q- 
cycle occupation density Cy(q) involves q distinguishable particles and n — q bosons 
and therefore the particle-hole symmetry is broken. 

• When g p ond > 0, g p ong starts below ^ ond since its slope at p cr is equal to 1 — 2p cr 
while g p ond has slope 1 — (3p cr and (3 > 2. Conversely, g p OQg finishes above g p ond 
since its slope at 1 — p cr is less than that of g p ond . 

Therefore in contrast to the Ideal Boson Gas, and the Mean-Field and Perturbed Mean-Field 
Models, as can be seen in the above diagram, the condensate density is not always equal to 
the long cycle density for the case of the Infinite- Range- Hopping Bose-Hubbard Model with 
Hard-Cores. However it remains true that the absence of condensation implies only finitely 
long cycles and visa-versa. 

In the following chapter we will consider again the Infinite-Range-Hopping Bose-Hubbard 
Model but this time without hard-cores to see if a similar result holds. However before doing 
so it is of interest to briefly consider the notion of Off-Diagonal Long- Range Order and its 
connection to Bose-Einstein condensation using cycle statistics. 





ong 



if and only if g p ( 



cond 



0. 
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Figure 4.2: ^ ond and g^ ong for /3 > 2 

4.7 Off-Diagonal Long-Range Order 

The one-body reduced density matrix for x, x' G A v may be denned as 

D njV (x,x') := {a* x a x ,) Ht c v = ]F^ trace ^>, + 
where for <p G H v , K xx i<f> = (e x /|^»)e x . 
Penrose showed that for x ^ x', 

lim D njV (x,x') = g p cond , 

n, V — ¥oo 

n/v=p 

that is, whenever Bose-Einstein condensation occurs, there is Off- diagonal long-range order 
as defined by Yang [97]. It has been argued and proved in some cases (see for example [81] 
and [34]) that in the expansion of D n v (x,x') in terms of permutation cycles, only infinite 
cycles contribute to long-range order. Here we are able to show this explicitly. 
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Using Theorem 2.2 from Chapter 2, we have the following 

n 

D njV (x,x') = y^c"(g; K x>x >) 

9=1 

where 

c n v {q- K x , x ,) = 1 trace W hc \(K XiX > ® I <g> I <g> . . . <g> I)U q e~^] . (4.18) 
Z v [n)v v> n > v L 

Note that this is equivalent to the expansion of cr p (x) in equations (2.14) and (2.16) from 
[81]. 

Applying the same argument as in Sections 4.3 and 4.4, one may show that 



lim c n v (q; K XtX .) = lim c^(q;K x ^) (4.19) 

n,V— ► oc n,V — >-oo 

n/v=p n/v=p 

where we take 



c£(g; K x , x .) = trace W hc {K x , x , ® I ® I <g> . . . <g> J)?7 g e 

Z v (n)V n i,n,v [ 

The only difference is that instead of equation (4.12), we obtain 



E 



E E E 

(K,C; • • • ;r 2 ,^ r 2 2 ;r 1 ,z r 1 i ](i°)|(^ <g> J® • • • ® /)t/ g i ) (4.20) 

whose treatment is similar but slightly more complicated, as detailed below. 

Let q > 1 and consider the case {ri, r 2 , . . . , r m } 7^ {1,2,..., q}. When 1 ^ {ri, r 2 , . . . , r m } 
we obtain inner products of the form: 

(il\K X)X >il)(j 2 ,h, ■ ■ ■ , j q \il,i°, ■ ■ ■ ,i° q ,ii) 
where jk 7^ for all k by the hard-core condition, implying the second term is zero as 
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j q ^ i\. On the other hand, when 1 G {ri, T2, ■ ■ ■ , f m }, then there exists at least one 
/ ^ {7*1, r 2 , . . . , r m }, yielding an inner product of the form 

(ji\Kx,x'il)(j2, ■ ■ -,ji-i,ii,ji+i, ■ ■ ■ ,j q \il,il, ■ ■ -,iq,ii) 

which also results in the second term being zero as = 0. Note that the above cases 

do not occur for q — 1. 

For the case {ri, r 2 , . . . , r m } = {1, . . . , q}, as before, the remaining sum may be bounded by 
a similar expression whose summations have slightly relaxed restrictions. Also the left hand 
side of the inner product is independent of i°, so again denoting s° = (1, 2, 3, . . . , q), we have 

V V V 

(4-20) < E E •• E E<[ r ™'C;---;^^;n,iy(s )|(^iS),i§,---,*-;,*-S> 



7 -l _i ,-2 _i 



,-1 U=Oi 



where [r m , ; ... ; r<2 , i r ; ri, i r ] (s ) 
has distinct indices 

and as there is only one possible value for each i\, i°, z'2, . . . , i° q giving a non-zero summand, 
we can bound above by 

< E E ■■• E £«J^S> = f- 1 E E«J^.^5) = v— 

,1 _i ,2 _i ;m =1 ,0_i ;k —1 ,0_i 

where k G [l,m] is the smallest number such that = 1, and for any x, x' G A v . Thus the 
entire sum (4.20) is bounded above by V m ~ Y . Therefore one can conclude the argument of 
Subsection 4.4, proving (4.19). 

Moreover, following the reasoning in Subsection 4.5, we can then check that for q > 1 and 
x ^ x', Cy(q; K XjX r) = 0, since for q = 1, (ei\K XiX >ei) = 0, and for q > 1 

{{K XjX , ® / <S) ■ ■ ■ <8 I)U q i\i) = {{K XtX ,e i2 ) (8 e h <8 • • • <8 e, <g> | e ix (8) e h <8 • • • <8 e, ) = 
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as the z/'s are all distinct. So we have that 



lim Cy(q\ K XjX >) = 

n, V— >oo 

n/v=p 



and that 

lim lim — V] cl(q;K XjX ,)= lim D n>v (x, x') = 

Q — >00 — foo \/ ' * 71, V— >oo 

n/v=p <?=Q+1 n/v=p 



thus proving the long-range order coincides with the condensate density. 
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Chapter 5 

The Infinite-Range-Hopping 
Bose-Hubbard Boson Model 

Summary 

In this chapter we study the relation between long cycles and Bose-Condensation in the 
Infinite-Range-Hopping Bose-Hubbard Model (without hard-cores) [14] and calculate the dens- 
ity of particles on long cycles in the thermodynamic limit. The argument is as follows, 
first we shall obtain an expression for the density of particles on cycles of length q in the 
grand- canonical ensemble, isolating q unsymmetrised particles. Secondly we prove that in 
the thermodynamic limit we may neglect the hopping of these q particles in the cycle density 
expression. Third we attempt to justify the application of the Approximating Hamiltonian 
Method upon the cycle density, proving results in the absence of condensation, and with a 
gauge-breaking term, in its presence also. Finally we apply simple numerical methods to the 
resulting cycle density expression to indicate our conclusion. We find, as in the hard-core 
case, that the absence of long cycles implies the absence of condensation and visa versa. 
In the presence of condensation, we argue that while the occurrence of Bose-Einstein con- 
densation coincides with the existence of long cycles, their corresponding densities are not 
necessarily equal. 
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5.1 The Model and Previously Derived Results 

The Infinite- Range- Hopping Bose-Hubbard model (without hard-cores) is given by the Hamilto- 
nian 

v v 

H ™ = 7v E {a: ~ < )(a * " a y ] + x E - ^ ^ 



2V 

x,y=l x=l 

where the magnitude of repulsion is controlled by the parameter A > 0. The zero-temperature 
properties of this model have been analysed for instance in [43; 44; 77; 40], mostly by applying 
perturbation theory to the ground state wave function or via computational methods. 

To observe the low temperature properties of this model, Bru and Dorlas [17] applied the 
"Approximating Hamiltonian" technique developed by Bogoliubov [12; 13] (see also [9; 10; 11] 
and [7; 8]) to this model. In this method one performs the following substitution for the 
Laplacian term of the Hamiltonian: 



1 V V 

y a * a v ~^ 5Z( fa * + ra *^ ~ y l r ' 2 

x,y=l x=l 



(some r G C) to obtain the approximating Hamiltonian: 

v v v 



x=l x=l x=l 



Introduce a gauge breaking source v G C to both Hamiltonians (5.1) and (5.2), by setting 
H v (u) := H v - ELiO^x + va* x ) and Hf F (r, v) := #£ pp (r) - ELi(^ + va D- Then for 
all n (with some constraints on u), one finds that the pressures for these Hamiltonians are 
equivalent in the thermodynamic limit, i.e. for large V one obtains the estimate 

< Py[Hy{v)\ - supp v [H^ P (r, u)] < 0(V^ 2 ) 

where p v [-] '■= p v [-}((3,fi) is the grand-canonical pressure. Henceforth the f3 and \i depend- 
encies are assumed unless explicitly given. 

The exactness of this approximation is discussed in [45; 61]. Recent applications of this 
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approximation to a continuum Bose gas model appears in [18] and [99]. This model in the 
presence of an additional random potential is considered in [30]. With this approximation 
technique, Bru and Dorlas managed to obtain the limiting pressure and showed that in some 
regimes Bose-Einstein condensation occurs. They proved the following result: 

Theorem 5.1 The pressure in the thermodynamic limit for the Infinite-Range-Hopping 
Bose-Hubbard Model, p(f3,fi) := limy^oo p v [H v ] , is given by 



p(/3,/i) = sup j - r 2 + - In trace j- + (c)e" 

r>0 [ P 



(5.3) 



where 

h(r) := (1 — fi)n + Xn(n — 1) — r(a + a*) 

is a single site Hamiltonian with creation and annihilation operators a* and a, and with 
number operator n = a* a. Note that now the supremum is attained over the set of non- 
negative real numbers. 

The Euler-Lagrange equation for the variational principle is 



tracejr ( o(a + a*)e /3/l(r) 



Moreover the density of the condensate is exactly given by 



! v 

? cond ^^^H ( a >v) Hv = T l 
x,y=l 



where r M is the largest solution of (5.4). 

Equation (5.4) can have at most two solutions. Clearly r = is always a solution. When /3 
is large enough, for certain values of fi a second non-zero solution may appear (see Figure 
5.1). So 77, := unless a second solution r > exists, in which case r M := r. 



92 



CHAPTER 5. I.R.H. BOSE-HUBBARD MODEL 



1.5 - 



1.5 - 



0.5 - 




0.2 
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0.6 



0.8 



Figure 5.1: Comparison of 2r with (a + a*)h( r ) with (3 = 4, A = 5 for the cases fi = 1.5 
(condensation) and fi = 5 (no condensation) 

Bru and Dorlas then obtained the properties of this model by using some numerical tech- 
niques to find this maximal solution to the Euler-Lagrange equation and use it to evaluate 
the pressure using (5.3). 

As may be seen from Figure 5.2, for sufficiently large (3, there may exist several critical values 
of ji which correspond to intervals of r M = and r M > 0. 

In addition Bru and Dorlas showed that Theorem 5.1 holds in the presence of the gauge- 
symmetry breaking term. In that case, the corresponding Euler-Lagrange equation has a 
unique non-zero solution r^iy). 

5.2 Cycle Statistics and the g-particle Hop 

From Chapter 2, we have established, in the grand canonical ensemble, that the density of 
particles on cycles of infinite length is 




q=Q+l 
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3 I 1 1 1 1 1 1 1 r 
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5 10 15 20 25 30 35 40 

M 

Figure 5.2: Plot of r 2 {— g^ ond ) versus \i, for (3 = 4, A = 5 

where c v (q) is the density of particles on cycles of length q. We intend to evaluate this ex- 
pression for the Infinite- Range- Hopping Bose-Hubbard model and compare it to the (already 
known) condensate density ^ ond . 

In order to proceed, it is best to reaffirm the notation used in this chapter. As before, the 
Hilbert space for a single particle on a lattice of V sites is T-L v := C v and upon it we define 
the operator 

k v = I- Py 

where P v is the orthogonal projection onto the ground state g v . We can define the free 
Hamiltonian acting on the unsymmetrised Fock space Tijiy) by the second quantization of 
the operator h v : 

Hy = dT(hy). 
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For bosons we consider the symmetric subspace of ^(Hy), denoted F + (H V ). 

The operator which counts the number of particles at site labelled x, denoted n x , is defined 
by n x = dT(P x ) : = dT(\e x )(e x \). Then the total number operator is N v := Y^x=i n x = dT{I). 

Then the Bose-Hubbard Hamiltonian acting upon J : + (1-L v ) may be expressed as 

v 

H v = dT(h v ) + \ ^2n x (n x - 1). 

X=l 

Note H v is defined for all particles (not just bosons), i.e. H v may also act upon J r ('Hy). 

Recall we fixed as the identity operator on "Hy\ and I as that for jF + (T-L v ). By Theorem 
2.1, we have the following expression for the g-cycle density: 

c v (q) = -L trace Hq ^ v [{U q ® I) e -«^-^] (5.5) 

where U q , v := U^®F + (H V ). 

Note that partitioning the partition function into its cycle structure is only possible if (2.2) 
holds, i.e. 

trace ^ [U n e^ Hv ] > 0, 

however here by the random walk formulation one can see that the kernel of e~ l3Hv is positive 
and therefore the left-hand side of this expression is positive. 

As is evident from equation (5.5), the cycle statistics technique of Chapter 2 results in 
splitting the symmetric Fock space J-+(H V ) into the tensor product of two spaces, an un- 
symmetrised g-particle space T-L^ and a symmetrised Fock space F + {H V ). Similarly to the 
previous chapter, it is convenient to have a notation to split the action of an operator on 
jF + {^-L v ) to those upon Hy and J Z + {H V ) individually. Write 

A (g) := A {q) g, J and dr'(A) := I (q) ® dV{A) 
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for any operator A on T-L v . In this fashion, the number operators applied to Hq sV are defined 
as 

N x := ® I and n x := I {q) <g> dT{P x ). 

Then in this notation we may define Hy on T-L qi y by 

v 

H$> = + dV'{h v ) + A + N x )(n x + N x - 1). 

x=l 

As in the hard-core case, we intend to drop the contribution of the hopping of the q unsym- 
metrised particles and see if this has any influence upon the cycle-density expression in the 
thermodynamic limit. To do so define the modified Hamiltonian by the following: 

v 

: = Jfe) + dT\h v ) + A J2( n * + + N x - 1) 

x=l 

so that Hy ^ = Hy — Py \ and then define the modified cycle density, that which neglects 
the hopping of the q distinguishable particles, by 



1 1 

c£(g) = ^r— trace n ^ 

^y V 



We wish to prove for each q that the cycle density c v (q) and the modified cycle density Cy(q) 
are equivalent in the thermodynamic limit, i.e. 



Theorem 5.2 



lim c v {q) = lim c^(q) 

V— >co V— >oo 



which implies that in the thermodynamic limit, we are able to disregard the hopping of the 
g-unsymmetrised particles in the cycle density. 

Proof: The proof of this Theorem is structurally similar to that of Lemma 4.1 from the 
previous chapter when we considered the Hard-Core model. However there are sufficient 
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differences between the two approaches to justify giving the complete proof. We wish to 
obtain an upper bound for the following modulus: 



trace % 



U ge -l3(Hv-i*N v ) 



trace % 



q.V 



Before proceeding, we first shall introduce a basis of the Hilbert space % q , v . Let {<pk}kLo be 
an orthonormal basis for jF + (7i v ). 

Let Ay^ be the set of ordered g-tuples of (not necessarily distinct) indices of A v and for 
i := (ii,i 2) ...,«,) G Av } let 



|i) = \k,i2, ■■-,i q ) ■= e h <g) e i2 <g) • ■ ■ ® e i? . 



Then {|i) |i G A^} is an orthonormal basis for 



(9) 



A basis for may therefore be formed by taking the tensor product of the bases of TLy^ 
and T + (H V ), so the set {|i) ® fe | k = 1, 2, . . . ; i e At 9) } is an orthonormal basis for 7^,, v -. 
For brevity we shall write 

|i;fc) := |i) ®<pk- (5.6) 



For simplicity, denote P := P^ and H := — /iA r v . We expand 



trace Hq V [U q e 



-f}{H v -»N v ) 



] = trace Hq V [U q e 



in a Dyson series in powers of P. If m > 1, the m term of this series is 



/■l /"SI PSm-1 

X m := /3 m / dsi / <is 2 • • • / ds m trace % 
Jo Jo Jo 



e _/3H(l- Sl ) p e -l3H( Sl -s 2 ) p . . . 



(5.7) 



Let P r := J <8> • 



r th placc 



J, so that P = YZ=t Pv Then 



97 



CHAPTER 5. I.R.H. BOSE-HUBBARD MODEL 



q q 



Xm = P m J2J2---Yl X m(ri,r 2 , ...,r m ) 



ri =l 7-2=1 r m =l 



where 



rl rsi rs m -i 

X m (ri,r 2 , . . .,r m ) = dsi ds 2 . . . / ds m trace H 
Jo Jo Jo 



q,V 



-PH{\-sx) p p -/3tf(si-s 2 ) p ... 



In terms of (5.6), the basis of %„ v , we may write 



X m (ri,r 2 , ■ ■ ■ 



"1 rs\ 

r m ) = Idsi ds 2 



JO 



E E-E 



fcO ,fc 



m ;0 im 



(5.1 



• ■ • (i™" 1 ; k m - 1 \e-f 3H{Sm - 1 - Sm) P r Ji m ; k m )(i m ; k m \e-^ Hsm U q \i°; k°) (5.9) 

where it is understood that the i summations are over Ay\ the set of ordered g-tuples (not 
necessarily distinct) of Ay, and the k summations are over the bases for J r + ('H V ). 

Notice that we may express 



e-^|i;fc) = e 



-Pq(l-»)s \- -PH l s 



\k) 



where 



H 1 ■= H v - fiN v + A Y,(K(K - 1) + 2AT> X ) 



x=l 



and N* := E^ =1 S x ,ij counts the number of particles at site x which are in g-space. Also, for 



any fixed r: 

1 V 

...,i r ,j,...,i q ;k) = y^2\[r,j](i);k) 

3=1 3=1 

using the same index replacement notation [ • , • ] as in Chapter 4. 
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Using these facts, a single inner product term of (5.9) may be expressed as 
(i; k\e~P Ha P r \y, k 1 ) = e"^ 1 "^ (k\e^ His \k') (i|P r |j) 



V 






-H)s 


V 




e -/8?(l- 




V 



(k\e-e Hls \k')J2(i\ju---Lm,...J q ) 

m=l 
V 



l qi3q 

m=l 



e 



\k } Si 1 j 1 . . . Si r j r . . . 5i q j q . 



Now if we sum over j 



Y$ k\e-? Hs P r \i; k')(y, k'\ = (k\e^ Hls \k') . . . , t r J r , . . . , i q - k'\ 

j >=i 

= y <* I e~^ s \k')J2( [r, Jr] (i); k' \ . 

j r =l 



Performing two summations for fixed r\ and r 2 we get: 

^^(i^^le-^Ji 1 ^ 1 ) (i^^le-^P,,!! 2 ^ 2 ) (i 2 ;k' 



r 



l l 1 =l i2 



e "^f ( ' +f) £ ^V^l^) (^|e-^- ,(i0) ^ 2 ) (h,^;^,^]^);^! 



7 1 =1 2 2 =1 



Thus (5.8) looks like 
X m (ri,r 2 , . . . ,r m ) 



e -/3g(l-At) /-i /-si ' 



. . . y E E E • • ■ E (^le-^-)^ i^> 

■'0 >/0 >/0 ,n :n .1 _i _i Am — i 



- - - <fc™|e-^— ^ E — ^ , - ra -^""^- 1C * 0> |fc ><[r- TO , «^*; - - - ; r- a , ; r-x, C*°) I ^A°> 
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— jdsx I ds 2 ... I ds m 



JO JO 

V V 



EEE-E<r 

:0 ,1 —i ,2 —i im = i 

1 V 1 - 1 'r 2 - 1 t r m 1 



r 2 ,e;r 1 ,^ 1 ](i°)|f/ ? i°) 



trace 



e _ /3 (l_ Sl )Hi° e - / 3( Sl -s 2 )H [ri ' I '-l] (i > e _ / g Smi? ['-m,«^;...;r 2 ,v 2; , 1 ,v 1 ](i ) 



Using the Holder inequality, for non-negative trace class operators Ai, A 2 , . . . , A m+ i we have 
that 

771+ 1 

[trace (AiA 2 . . . A m+ i) | < trace |Ap4 2 .. < 1 [ (trace A^ ■ 



j_ 



k=l 



where Yl=i h = h Pi>®- 



Set pi = T ^ 7 , p 2 = 



j • • • j Vm 



— — ,p m +i — —■ Taking the modulus of the above 



trace 



trace r+{u v ) 



e _ J 8fli°(l- ai ) e -j9ff lri,<J i 1(,,}) (*i-*3) e _ /3H [-.^„;-;''2,^ 2 ;n,4 1 ](i ) (Sm) 



< trace 



trace t + ( Hv ) 



trace T+{nv) 

-/3H l m ' Tm ' "' '''a"rai'"i'*ri 



r 1 i^CKl s l~ s 2 
Q H -[ r l>«f 1 K» ) 



But since the trace is independent of the V — q sites {i°, [r 1; [r 2 , ^ ; r 1; ^J(i°) 

. . . , [r m , \ . . . ; r 2 , ijL; ri, 4j(i )}, the product of all the trace terms above is equal to 



trace t+(h v ) 



-/3H l 



with 1 = {V - q + 1, V - q + 2, . . . , V}. 



This is independent of the i° and i summations, so we need only consider 



v v 



EE-EE* 



. ;r 2 ,i^ -n.ii ](i°)|(7,i°). 



l™ =1 iO 
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Fix the values of il , . . . , i™ . We intend to show that 

i 1 1 i 2 1 1 1 m 

^([r m ,z™;...;r 2 ,^;r 1 ,^ i ](i°)|f/ (? i ) = l. 
i° 

If {ri,r 2 , . . . ,r m } ^ {1,2, . . . ,q}, then |[r m ,z™ ; . . . ; r 2 , i 2 n ; n, ^J(i°)> is of the form 

J2l ■ ■ ■ 1 jni ) + • ' • 5 ^n2' J". 2 + l? • • • ) jVt 3 ; ^ri3 + l) • • • 5 ^714) jVl4+l) ) 

where {ni,n 2 , . . . } is a non-empty ordered set of distinct integers between and q. This 
vector is clearly orthogonal to U q i° except for the single choice of 

i 1.72 j • • • j jni— 1) jnii ■ ■ ■ i 3mi Jn2+1) ■ ■ ■ > jna—li jn%i ■ ■ ■ i 3nzi jm+li i jl) ■ 

For the case {n, r 2 , . . . , r m } = {1,2,...,?} notice that |[r m ,i™; ... ; r 2 , z^; ri, «JJ(i )) is 
independent of i° so we may take it to be 

|[rm,C; ••• ;r2,i£,;r-i,iy(s )) 

where s° = (1, 2, 3, . . . , q). For each choice of ,if ,. . . , i™ there exists only one possible 
i° E A ( v ] such that 

([r m ,i£; ... ;r 2 ,z^ ri , z i i ](s°)|^i )^0. 

So we may conclude that 

EEE-E ([r m ,C; ... ;r 2) 4;r 1) <J 1 ](i )|Z7 ff i > = V" 

:0 ,'1 —1 ,;2 —i jm =1 

and by applying this, we see that the modulus of (5.8) may bounded above by 
\X rn (r q ,r 2 , ...,r m )\< trace F+{Uy) [e"^'] e^ 1 -")^— 



m\V m 

EEE-E ([r m ,C m ; ... jr^^r^^KOlC//) 



:0 ,1 —1 ,-2 —1 =1 



trace 



e -/J?(l-M) J__ 

m! 
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which is independent of 7*1, . . . , r m . Hence the modulus of (5.7), the m th term of the Dyson 
series, may be bounded above by 



\X m \ < /3 m ^--- ^ \ x m(r q ,r 2 ,...,r m )\ < trace 



ri=l r m =l 



-PH 1 



n mor 



ml 



Noting that the zeroth term of the Dyson series is 



Xq = trace ^ v [U q e ^ ] = trace 



H 



q,V 



Uqe -P(H^-^Ny) 



we may re-sum the series to obtain 



trace Hq>v [U q e~^-^] -trace^ [u^v-^v 



< trace t+(h v ) 



trace f+{h v ) 



-PH 1 



-PH 1 



q m fi r 



e^(l-e"^). 



//?! 



Thus 



K(q)-c$(q)\ 



V 



trace « v [C/,e" 



■P{H v -nN v ) 



] - trace w [Z7,e _ 



/3(^ ) - A1 7V v )l 



< 



e^3M(i _ trace jr +(Wv) [e"^ 1 ] 



I" 



Finally considering the ratio of the traces, since H v — fiN v — H 1 = A $^x=i(-^x(-^x ~ 1) + 
2N]n x ) > 0, the second fraction is not greater than 1, implying 



<(?)-c£(?)|< 



e^(l -e - ^) 
V 7 



which goes to zero in the limit V — > 00, as desired. 



□ 
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Corollary 5.1 The density of cycles of finite length in the IRH Bose-Hubbard model may 
be expressed as 

a = Y^ e -p iq - fl)q Hm trace F+mv) exp { - /3(2Agn! + H v - fiN v )} 

short ^ ^oo trace T+{Hv) exp {-(3(H V - fiN v ) } 

where ri\ is an operator which counts the number of bosons on the site labelled 1 . 

Proof: Note that we may easily re-express the modified cycle density expression: 

e -/8(i-A*)9 ~ 
c£(q) = ^ v 2^ • • • 2^ 2^{h, i2, ■ ■ ■ , i q ; k\ 

«i=l ig=l k=l 



V V 

n x 

x=l x=l 

U q \ix,i 2 , ...,i q ;k) 



exp I -P dT\h v ) + A J2(n x + N x ){n x + N x - 1) - ^ 



q | "1 j "2i • • • i "qi 

V V oo 



= ^ v 2^ • • • 2^ Zj(h, k, ■ ■ ■ , iq\U q \ix, z 2 , • • • , i ff > 

^ H=l ig=l fc = l 

{/ V q q V 

-f3 I cZT(/i v ) + A^(n x • .,•)(», + ^ - 1) - H ^ ) }> 

\ X = l j = l j = l X = l 

and as 12, ■ ■ ■ , % q \%2i «3, - - - , i q , h) 7^ if and only if i\ — %i — ■ ■ ■ — i q :— i then 



-0(1-/09 



exp \ dT(h v ) + \^2(n x + qS ix )(n x + g<5 ix - 1) - // ^ j J \k) 



e -j9(l-A*)« 



^— - ^ trace jr +( „ v) 
^ v i=i 

exp < - dT(h v ) + \{rii + + g - 1) + A^n x (n x . - 1) -jU^ra x 

^ \ x — l X = l / 
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e -P(l-H)q 



(3 dT(h v ) + X(m + q){n l + q- 1) 



trace t+{h v ) ex P 




x=2 x=l 




since the trace is independent of the basis chosen. Hence we obtain 



c${q) = e 



trace t+(h v ) ex P { _ P{^M n i + H v - fiN v )} 
trace t+{u v ) exp { - (3(H V - fJ,N v )} 



where n\ is an operator which counts the number of bosons on the site labelled 1 of the 
lattice. The result follows by taking the thermodynamic limit and summing over all q. □ 

5.3 The Approximating Hamiltonian Method 

In Corollary 5.1 we derive the following expression for the density of cycles of finite length 
in the IRH Bose-Hubbard model: 



where n\ is an operator which counts the number of bosons on the site labelled 1. Unfor- 
tunately we are unable to treat this expression analytically. Instead we apply the Approx- 
imating Hamiltonian substitution to this expression and conjecture that the approaches are 
equivalent in the limit. 

If one were to substitute H v for the approximating Hamiltonian Hy PP (r^) (where again r M 
is the maximal solution to the Euler-Lagrange equation (5.4)) into the right hand side of this 
expression, one would obtain: 




q=l 



oo 



lim 

V— >oo 



trace t+^w) exp { - /3(2Agni + H v - fiN v )} 
trace t+(h v ) ex P { - Pi H v ~ ^N v )} 




(5.10) 
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where 

h q (r) := (1 - fi)(n + q) + A(n + q)(n + q — 1) — r(a + a*) 

is another single-site Hamiltonian (note that ho(r) = h(r)). We conjecture that (5.10) gives 
the correct expression for g^ hovt . 

Moreover the fact that a state corresponding to H v in the thermodynamic limit may be 
shown to be a convex combination of one-site product states of the form 

trace jr r C ) z~ phq (r " ^ A 
v ' trace j- +(c) e-^o( r M) 

supports this conjecture. We shall prove this conjecture for those values of /x such that 
r M = 0, but unfortunately are unable to do so when r M > 0. However we can prove an 
slightly weaker result with the addition of a gauge-breaking term. 

Let Cv{q, v) be the density of particles on cycles of length q for the gauge-symmetry broken 
Hamiltonian H v {y) and set c M (g) := limy^oo c^(q). 



Theorem 5.3 For ji such that = 0, we have 



trace jr + ( C )e 0h "^ 



trace jr +(C )e-^o(°) 



More generally for any /j, 6 R, for a fixed v > there exists a sequence u v — > v as V — > oo, 
independent of q such that 

trace ^ +(qe -^ 9 (r M H)-K<M-a*)] 



Proof: For convenience, denote H q ^ v := 2\qn\ + (q — jj)q + H v — fiN v . Due to Theorem 5.2 
we need to consider 

trace t + {h v )[^ H "' v ] 



c fl (q) = lim 



v^oo trace t+(h v ) [e~P Ho - v ] 
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Recall that 

1 V V 
Hv = zv Y ~ a 7>( a * ~ a v"> + A Y n ^ Ux ~ ^ 

x,y=l x=l 

Motivated by the occurrence of the site-specific operator ri\ in the numerator of the expression 
of Cy(q), we expand the expression for H v — fiN v to isolate the operators which apply to the 
site labelled 1: 

H v - fiNy = (1 - n)ni +\n 1 (n 1 -l) + y-y^a*-^^2a x + H v 

x^l x^l 

where 

„ i 

Hy := (1 - fl) ^ n x ~ y Y + X Y ~ ^ 

xj^l x,y^l x^l 

is an Infinite- Range- Hopping Bose-Hubbard Hamiltonian for A v \ {!}. By denoting 



h q , v := (1 - n)(m + q) + A(m + q)(m + q - 1) + f 

we may then write 

H q ,v = h qtV - yY a * x ~ yY ax + 

x^l x^l 

Note that h qv — > h q on J>(C) as V — > oo. 

We intend to completely segregate the Hamiltonian H q ^ v into two individual parts, one which 
operates solely upon the site labelled 1, and the other which applies only to the remaining 
V — 1 sites. What prevents us from doing this immediately is of course the "cross-term" 

x^l x^l 

Motivated by the Approximating Hamiltonian technique, we shall substitute this term with 



106 



CHAPTER 5. I.R.H. BOSE-HUBBARD MODEL 



for a certain co-number R. Without loss of generality, we may take R to be a non-negative 
real number. Fixing 

h q , v {R) ■= hq jV - R{a x + a*) 
then the resulting newly approximated Hamiltonian may be expressed as 



In the arguments that follow, we shall either take R to equal r M in the variational principle, 
or a variable depending on V which tends to r M in the limit. 

5.3.1 Case 1: Values of /i such that r M = — the absence of condensation 

First we shall state and prove the following: 
Proposition 5.1 For all ji Gi such that r M = 0, 



H™(R) = h q , v (R) + H, 



lim 

V— >oo 



trace T+{Hv) exp { - PH* p v p (0)} 
trace t + (h v ) exp { - f3H q , v ] 



1. 



Proof: Using the Bogoliubov inequality: 



(A-B) B < In trace e A - In trace e B < (A- B) A 



(5.11) 



for any R we obtain 




< In trace exp { — (3H qtV } — In trace exp { — f3H^(R)} 




+ h.c. (5.12) 
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Taking the left-hand side, since H^(R) is a sum of two Hamiltonians which act upon 
different Hilbert spaces, traces and therefore expectations may be easily de-coupled, so one 
may see that 

which is zero since H v is a gauge-invariant Hamiltonian: ^Ez^i a*^j _ = 0. 

On the right-hand side, note that {a\)H qV = 0, also due to gauge invariance. Therefore we 
may simplify (5.12) to obtain 

-PR{ ai + al) hqv{R) 

< In trace exp { — (3H gtV } — In trace exp { — f3H q ^ p (R)^ 

*'M^)L +h * (5 - i3) 

For this case we shall take R = r M = 0. We therefore obtain 
< In trace exp { - (3H q , v } - In trace exp { - (3H^ V F (0)} < (3 (a 

Now by the Schwarz inequality 



mi,"" 

(5.14) 



E 



x£l a x 
V 



II 



a 1 a 1 



H, 



q,v 



V 2 



q,V 



V 



To consider this let H s q v := H qtV + s '^ a * x ] p ax and H* v := H s qv -2Xqn l . Using the Bogoliubov 
inequality (5.11), with A = —/3H* V and B = —/3H* V , so that A — B = 2j3\qn\, one has 

2(3\q(n 1 ) H s v < lntracee _/3 ^V - In trace e~^V < 2p\q(n 1 ) s . v = 2f3\q — — (5.15) 
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The last equality may be seen due to the fact that the system H^ v is invariant under 
permutation of the sites of the lattice. This identity implies the following: 



i) With s = we get (n 1 ) Hq v < (N) Hv /V, thus (n 1 )u q v is bounded and in the limit 
-»■ 0. 



V 

(ii) 



< — In trace e _/ "V In trace e _/3i V < 



v v - y y 

indicating that in the limit, the pressures are the same for the two Hamiltonians. 
By using Griffith's Lemma we may see that the condensate densities (the derivat- 
ives with respect to s at zero) are both equal to zero (since we are considering the 
case r M = here). That is 

Urn Pip' ) = Um ( E tP' ) ==0. 

v -°°\ V /ff„v V -""\ V I Hv 

Using these facts, one may see that the right-hand side of (5.14) goes to zero in the limit 
and we can conclude that 

trace t+{H v ) ex P { ~ P H q,v) 
lim - ( — t^-v = 1. 

k->oc trace T+{Uv ) exp { - f3H^/(0)} 

□ 



Using Proposition 5.1, one immediately obtains the desired result: 

at \ r at \ 1- traCe ^+(W V )eXp { - PH qtV } 

(?{q) = hm c^(q) = hm p ^—!r 

v^oo V ^oo trace jr +(nv) exp { - /3# ,v j 

. trace j- (w ) exp { - f} Hi^ (0) } 

= hm ? — — v 

v->oc trace T+{Hv) exp { - (3H^ V P {0)} 

lim trace:F +(^- 1 ) exp { ~ ^-^M trace HQ ex P { _ 

trace exp { - /?i?V} trace j- (c) exp { - f3h^ v ) 

trace jr(c) exp { - /3h q (Q)} 
trace j-(c) exp { - /3h (0)} 

in the limit since h qtV = h q , v (0) — > h q (Q) on .F((C). 
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5.3.2 Case 2: For any /i 

Consider the case when > If we insert R = r M in the constraint inequality (5.13) then 
its left-hand term is strictly negative, but its right-hand term is strictly positive, rendering 
the previous argument useless here. 

We therefore introduce a gauge-breaking term v Y2 X a x + v ^2 X a* into the Hamiltonians H q v 
and H^{r^). Without loss of generality we may assume v to be real and positive, so denote 

v 

H q , v {v) := H g!V - v ^2(a x + a* x ) 

x=l 

and its corresponding approximation as 

v 

H™{R, u) := H^(R) - u + <)■ 

x=l 

Again we wish to separate this Hamiltonian into parts, one acting upon the site labelled 1, 
the other on the remaining V — 1 sites. If we define: 

h q , v (r, v) := h q>v (r) - v(a x + a*) 

= (1 - + q) + A(m + q)(n 1 + q-l) + f-(r + i/)(ai + a{). 

and 

x^l 

then we may write H q ^ p (R, v) = h q>v {R, v) + H v {v). Denote liniy^oo h q>v {r, v) = h q (r, v) on 
J-" + (C), i.e. obtain a "gauge-symmetry broken" single site Hamiltonian 

h q (r, v) := (1 — pi)(n + q) + X(n + q)(n + q — 1) — (r + u)(a + a*). 

We shall first prove the following proposition: 
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Proposition 5.2 For each v > 0, there exists a sequence \y v E K : v v G [u, v + 1/y/V]} 
independent of q such that 

trace jr {H ) exp{- (3 H q:V (i^ v )} 
hm — — — f — - — -Y = 1 

v^oo trace t+{h v ) ex P {-P H qyM v v), v v )\ 

where r M (i/) > is the non-zero solution of 

2r = (a + a*) hq ( rjV ). (5.16) 

Note that lim^o r niy) = r 7t> the maximal solution to (5.4), i.e. the positive square root of 
the condensate density. 

Proof: There is no immediate correlation between the chosen R and as yet. For each 
v > take a sequence v v which tends to v as V — > oo. Then using the Bogoliubov inequality 
again, we obtain 

pL(^£*-r)) +h.c. 

< In trace exp {— f3H q ^ v (u v )} — In trace exp {—(5H^(R, z>y)} 

<pL(p^-R)) +h.c. (5.17) 



H q,v(yv) 



As above, the left-hand side may be reduced to 



Ei^l °*x D \ \ _ / \ / ^s^l °*x 



°1 77 Rj) = ( a l)h q y(R,u v ) ( 77 -~R 

V v / / #£pp (jw) V 

If we replace R with the term 

\ / ff V (^) 

then the left-most side of (5.17) is identically zero and we get that 

< In trace exp {—j3H qtV (v v )} — In trace exp {— f3H^(r^ v (i / v ), v v )) ■ 
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Hence 



liminf trace exp{-£ff g , y K)} ^ 

v^oo trace exp {-(3H^/(r- v (v v ),vv)} ' 



(5.18) 



Now considering the right-hand side of (5.17) for any \i. Using the Schwarz inequality as 
before: 



v 



H q ,v(») 



v 



v 



*«wi%h(^-*))L, 



H q y{y) 

V 



(ni>fl,, v (^(Vo) ff(iHl/) ) + y 



(5.19) 



where we have taken 



5 :~- 



Mt a -- vR )- 



Again in order to consider this, insert H qiV (u) and H q>v (b>) := H q ^ v (v) — 2Xqrii into the 
Bogoliubov inequality, to obtain 



2f3\q(n 1 ) Hqv(u) < In trace v{v) -In trace e -/3 ^ v(i/) < 2^Xq(m) S y(u) = 2f3\q 



V 

(5.20) 



which implies the following facts: 



'i) {ni)H q>v {v) < {N)a qv[v) /V, and hence <ni) y (y) ->. o as V ->■ oo. 



(ii) Since H q ^ v (v) = (q - n)q + - //iV, 



r • 



< 



0(9 - A*)? + I i ntraC ee-^vM-^v)_I lntraoe e-^.vM < 2 ^ 



U V 



V 



V V 



with which one may show that in the limit, the pressures are the same for H v (v) 
and H q>v (u): 

lim p v [H v (u)] = lim Pv [H q . v {v)}. (5.21) 
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Using fact (i) from above, we find that the only term we need yet be concerned with on the 
right hand side of (5.17) is the first term of (5.19), whose behaviour in the thermodynamic 
limit is still unknown: 



To deal with this we shall take R to be the following: 



\ x=l 



H 9 ,v(") 



so that one has 

v i v 



x=l X x=l 1 H qiV (v) 

Now we shall state and use some lemmas, which are proved later: 

Lemma 5.1 For fixed v > 0, a positive integer q, and 8q is defined as (5.22), then there 
exists a sequence {y v G R : v v G [v, v+l/y/V}} independent of q which tends to v as V — >• oo, 
such that for large V we have the approximation 



for some constants u, w > and M, independent of v and q. 



Using this lemma, for large V and fixed q, we obtain the following estimate 

coiisti 

In trace exp {-/3H q>v (u v )} - In trace exp {-(3H^ v p (r+ v (v v ), u v )} < 

where u v — > v as V — > oo, implying that 

trace exp {-^H q v (v v )} 

hmsup p ' — -y < 1. (5.23) 

v^oo trace exp {-PH^(r^ v {y v ), v v )\ 
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Lemma 5.2 For a fixed v > and for any sequence {zv} which tends to v as V — > oo, then 



lim r^ v (v v ) = lim r+ (u v ) = r^v) 



where r^iy) is the unique non-zero solution to the Euler- Lagrange equation (5.16). 



For clarity, it is best to use the following short-hand to complete this argument: 



a v = trace r{Hv) exp { - (3H q:V (u v )} 
b v = trace HHv) exp { - /3H^(r^(u v ), v v )} 
c v = trace T{Hv) exp { - pH£?(r~ v (u v ), u v )} 
d v = trace T{Uv) exp { - fiH^ v p (r+ v (v v ), u v )}. 



The penultimate step is to prove the following: 



lim — = 1 and lim ^ = 1. (5-24) 

V->00 Cy F->00 dy 



Considering the first, note that 

by trace ^ (c) e-^M"vW) 
c v trace ^(Qe-^^Vv-^v),^)" 

Once again using the Bogoliubov inequality (5.11) with A = —f3h q (r~ v (v v ) ) v v ) and B = 
—f3h q (r fl (h'y), v v ), then A — B = f3{r^(u v ) — r~ v (vy))(a + a*) and we obtain 

/3(rvK) - r~ v (u v )) (a + a*) hq{ril(yv)>uv) 

< In by — In c v 

< P{r„{v v ) - r- v (vy)) (a + a*) Mr - y(lvW) 

and as V — > oo, both the left and right hand sides go to zero, implying the result. A similar 
procedure may be used to show the second. 
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To complete this proof, we wish to derive 

lim — - = 1. 

V->oo b v 

with the use of the following information derived above (see (5.18), (5.23) and (5.24)) 

lim inf — > 1 lim sup -j- < 1 

lim — = 1 lim = 1. 

V->00 Cy V^OO dy 

Taking the supremum limit, one may see that: 

| < limsup^^ 
%~ hm^oofc 

The infimum limit follows similarly, proving the proposition. □ 



,. °V^ r dv y dy- . 11111 bU P^CXD Ty 

hm — < hmsup — = hmsup -r- < r — — < 1. 

V^oo by V ^oo by y^oo ^ hmy^oo j" 



With the assistance of Proposition 5.2 we have our result: 

v-+oc v^oo trace t+{u v ) exp { - pH 0jV (u v )j 

trace t+Wv) ex P { ~ l 3H ^ P ( r ^A 1/ v)^ u v)} 
v^oo trace r+{H v ) exp { - pH^{r^ v {v v ), u v )} 

trace exp { - /3H v (u v )} trace j- (c) exp { - Ph qtV (r^ v (i/ v ),i/ v )} 

v ^°° trace t+(h v -i) ex P { _ AW} trace j- (c) exp { - ^ v (r^ v (v v ), u v )} 
trace jr(c) exp { - f3[h q (r ^(u)) - u(a + a*)}} 
trace j-(c) exp { - f3[h (r^(u)) - u(a + a*)]} 
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5.3.3 Proofs of Lemmas 
Proof: of Lemma 5.1 

Considering the term (8q5 ) h ^ v ^ in a similar fashion to Appendix 1 of Dorlas and Bru [17], 
note that since [5 , 5q] = 1 we may write 

( 6 o S o)h cuV (v) = 2 (KA)>// q , v H - 1 
where {X, Y} := XY + YX. Using the spectral decomposition of the Hamiltonian H q ^ v [y) : 

H q , v {v)il) n = E n tp n 
and denoting A mn = (i/; m , 5 ,i/; n ), one may get 



m,n y ' 

+ ie-^^Ml £ |A mn | 2 (e-^ + e"^) /?(£ m - #n) 



6 

< 2(5 AK,vM + f <[*o, [^,v(^),5 ]])^.vH- (5-25) 

where the Duhamel inner product ( • , • ) is defined as follows: 

i rP 

(A,B) H = — trace [A*e~^- s)H Be- sH ]ds (5.26) 
Jo 

and Z = trace e~^ H . 

Denote Co = ^= a x . Consider the second term of (5.25) while using the following identities 
(for some x — 1 . . . V) 

2 1 

[ C 0, [«x- Cq]] = ^(2^ + 1), [Co, K, Cq]] = — , [Cq, [CqCq, Cq]] = 1, 

[4 [ c o, c o]] = [ C 0> [ c o, c o]] = 
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one may evaluate 

[5*,[H q , v (v),5 ]] = [c* ,[H q , v (v),c \] 

= [cq, [2Agni + (q — fi)q - c* c + (1 - A - fi)n x + An 2 + uy/V c* Q + uW c , c }} 

v v 
= 2\q[c* , [m, Cq]] - [eg, [c c , Cq]] + (1 - A - fi) J^[cg, [n x , c ]] + A J^[cg, [n 2 , c ]] 



Z/V V[Cq, [Cq, C ]] + I/VF [cj, [C , Cq]] 

2AV + g 



-/i + A + 2A- 



V 



Some useful inequalities and bounds 

Note the following operator inequalities: 



2V ^ v x v ' y y ' ~ ^ V 

x,y x x,y 

a; ^ ' a; 

z/Vy(c + c*) < z/ 2 c*c + V < u 2 N v + V. 

From these it is clear that the Hamiltonian with sources, H q v (v), is superstable for fixed 
q > and A > 0, i.e. 

TV 2 

H q , v {y) > A-rr - (A + /i + z/ 2 )iV v - V. 

By the Bogoliubov inequality, since H v (f3, /i — 2gA, i>) — H qv ((3, fi, v) — q{q — ^) = 2?A^ i=2 nj, 
one has that: 

X>) <(E^) • (5-27) 

We can find an upper bound for the expectation of the number operator with respect to 
H q v {v) which is independent of q. 

KUH^U + lt^ +(E«* 



i=2 / H qiV {y) \ " / H v (i/) \ i=2 / H v (p >fJ ,-2q\,v) 
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using (5.20) and (5.27). Now 
Y_ \ 

i=2 1 H v (l3,[i-2q\,v) 

since (N V } H rg -\ is monotonically increasing in \i. Hence we have a bound independent of 
q: 

The pressure p v [H v (i/)] is convex in fj, G M and i/Ei If one fixes the real numbers [Iq and 
z/o > 0, then by superstability there exists a uniform bound M := M(/3,/i) > such that 
(using (a x ) Hqv{v) = (a* x ) Hqv{v) and (5.28)): 

(c.f. (5.15)) with ji < fi and v < u . 

Resuming the proof 

Then we have 

(<W,f/ 9 , v h = (<5o, <>o)fl-,,vM + 3 I -/^ + A + 2A ^ J - 1 



< (*o, *o)h 9 , vW + f (-/i + A + 2A (2M + | 



It remains to consider the inner product. For convenience denote f q ,v{v) := (<5o, 5o)fr 
If v G C, then H q v (v) = H q v — \/V(vc Q + i>Cq) and so we would have the following: 

9 d 

However as v G JR., we may use polar co-ordinates to write the following: 

(5.30) 
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We first wish to show that fq,v(y) is bounded by a constant independent of q. For fixed V 
and v, recall we have 



2VM d 2 tTT , 



from (5.32). So it remains to consider the latter term: 



2 



2 



I 3 l ^-iPv\ H iA u )\ = (( c o + c *o)*> c o + c o)H q , v{ ») ~ ( c o + Co/^vM 
< (c Q + Cq, c + %) Hq<v ^ 

and using the fact that (A, A) < \{A*A + AA*), we may obtain 



< (( c o + c* )(c + c* )) Hqv{u) 

= (c* c + c c* ) Hqv(i/) + 2\(c c Q ) Hqv{u) \ 

< ^(2M + 1 + y/M(M + l)) 



from (5.28). To summarize: 



/^)<^; + ^(4M+l). (5.31) 

Note however that this procedure does not yield a satisfactory bound for our purposes, since 
we desire a bound of order V a , where a < 1. Taking a different approach, fix q, consider 
(5.30), multiply both sides by v and integrate: 



l' + C 



for [v, v + e] C [0, u }. 
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By (5.29), we have for all q £ N, \i £ K. and v < uq that 



dv 



9 MH q ,vC)] = ^y( c o + c o) Hq , vH = ^ 



o) Hq , vH < 2VM (5.32) 



and so 



d 











< 


2{y + e)\[M 


+ 


2uy/M 


V 









2(2u + e)^M 



with which we obtain 



(5.33) 



From (5.31), \f q ,v(v)\ < CV where C is a constant independent of q and v' £ \y,v-\- e]. Now 
let e = l/^ 1 / 2 . Denote 

oo 

F Y (u) :=^e-%, v (^). 

Since f q>v {v) is uniformly convergent in z/ (see (5.31)) and each term is continuous, then 
F q (u) is continuous. From (5.33) may obtain: 



v'F v {y')dv' < 



(5.34) 



By the Mean- Value theorem, there exists some v v £ [v, v + V 1 ^ 2 } such that 

v'F v {i/)dv' = -LF V K) 
v 



(5.35) 



which in combination with (5.34) yields 



v v „ / x VM / 1 \ 1 VM / 1 



For any positive integer q, since e q f q>v {v) < F v (v), then 



e q UvM < — — I ^ 



/5 2# 
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so that 



fq,v( V v) < W~ 



/3 \ v v 2v y 

Then as V increases, we have a sequence {u v e R : i/ v G [i/, v + 1/y/V]} satisfying (5.35) 
which tends to v > as V — > oo, so that for large V we have the estimate: 



„ , , eVM / /77 1 



□ 



Proof: of Lemma 5.2 

Before proceeding, we need to show the following: for fixed v > 

lim p v [H v (v)\ = lim p v [H v (u)}, 

V— >oo V— >oo 

i.e. a single site's contribution is irrelevant in the thermodynamic limit. Recall that 



Hy = ^V ^2^~ a *v^ a ^~ a y^ + X ^2 n ^ nx ~^ + \ V(V- 1) ~ W *52 n 
The corresponding pressure may be expressed as 

p v [H v ] = In trace jr +(Hv) exp |-/3^ v | 
1 



r / i y_1 

— In trace t+{-h v -i) ex P \ ~ f 3 ( 7^7 ^ ~ a V>^ ~ a ^> + A n ^ n * ~~ ^ 

V \ T 11— 1 T — 1 



(^-) 2 P1 .. 1 [^_ 1 (/3(V).A(^),^ 



£,3/=l x=l 

v-i 



V(V-l) 



where abusing notation temporarily we have explicitly included the parameters of the IRH 
Bose-Hubbard Hamiltonian, i.e. we write the pressure of (5.1) as p v [H v ] = p v [H v ((3, A,//)]. 
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Then in the limit, with the use of the Bogoliubov inequality, one may verify that 



fy -1\ 2 

timp v [#v] = ^lirn ( — — J p v _L (^r) , A (yzr) ,/i - y^iy)] 



lim PvrffTv-] =p(P,n). 



Now proceeding to prove this lemma, recall that we chose 



r »A V ) y 



where H v {v) = H v — vYlx^i( a x + a x ) is the gauge-broken IRH Bose-Hubbard Hamiltonian 
on all sites of the lattice barring the site x — 1. Fixing a value of v > 0, there exists a unique 
r ^( z/ ) > as the solution to the Euler-Lagrange equation (5.16), i.e. 

2?>0) = (a + O/^Mi/),!/)- 

The pressure p v [H v (u)] is convex in v and its thermodynamic limit is differentiable for all 
v > 0. By Griffith's Lemma, we have 

lim iv[# v -(z/ v )] = ^- lim p v [£T v (i/)]. (5.36) 

The left hand side of this evaluates to 

d ~ 1 2 

As shown above, we have that 



lim p v [H v (u)] = lim p v [H v {v)\ 

V— s>oo V— i-oo 

:= -r^) 2 + - In trace exp {/?[(// - l)n - Ara(ra - 1) + (r^v) + i/)(a + a*)]}- (5.37) 
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so the right-hand side of (5.36) will become (also using (5.37)) 
d 

— lim py[Hy{v)] 
av v^c 



d_ 

dv 



-r^iv) 2 + — In trace exp (/3 [(/i — l)n — \n{n — 1) + (r^u) + u)(a + a*)] ) 



c/t* (^) 

^ ((a + - 2r^(v)) + (a + a*)h,( r# , MiI/) 

(a + a*)h,M")>") = 2r>(z/). 



as desired. 
Similarly, taking 



1 / V 



x=l I H q V (y) 



where H qtV {v) = H q ^ v — v Yl^=i( a x + a *x)- Label the corresponding pressure for this Hamilto- 
nian as p v [H q>v {v)]. Recall the expression (5.21) that we previously derived: 



lim p v [H q>v (u)) = lim p v [H v (u)}. 

V— >oo V— >oo 

As above by Griffith's Lemma, we have 

1 d 1 d 

lim r+(i/ v ) = - lim —p v [H q>v (v v )} = -— lim p v [H v (v)\ = r^v) 

V— >oo 2 V— >co Cll/y 2 aV V->oo 

as above. □ 



5.4 The Absence of Condensation 

Using the first statement in Theorem 5.3 we can rigorously show in the absence of condens- 
ation that all particles are on short cycles. 
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Lemma 5.3 In the absence of condensation, i.e. for those \i such that r M = 0, the density 
of particles on short cycles equals the total density of the system, that is: 



Proof: For those values of \i such that r M = 0, i.e. in the absence of condensation, first note 
that the density (from 5.3) may be expressed as: 

d d 1 

trace j-( C ) [nexp { — (3[(1 — fi)n + \n{n — 1)]}] 
trace jt(- C ) exp { — (3[{1 — n)n + Xn{n — 1)]} 

Similarly, from Theorem 5.3 one immediately obtains: 

trace jr (C ) exp { - (3[(1 - /i)(n + q) + \{n + q){n + q - 1)]} 



(f(q,0) 



trace exp { — — ji)n + \n{n — 1)]} 



Label the denominator S := trace jyc) exp {— — fi)n + Xn(n — 1))}. The operator n in 
this context counts the number of particles on the site, so in terms of a basis of occupation 
numbers, it has eigenvalues k — 0, 1, 2, Summing over this basis 



^ oo oo 

E cM (<7, 0) =- E ex P { " ^W 1 "/*)(* + «) + A(* + ?)(* + g - 1)] } 



g=l g=l fc=0 



and shifting the sum 

1 oo oo 

=i E E ex p i - w - ^ k + Xk{ * k - !)] } 

q=l k=q 
1 oo fc 

=s E E ex p i - 0K 1 - ») k + Xk ( k - x )] } 
1-1 /c=i g=i 

1 oo 



fc=i 
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=— trace j-( C ) [nexp { — — fj)n + \n{n — 1)]}] 

=p(h). 

Therefore the absence of condensation implies that the sum of all finitely long cycle densities 
equals the system density. □ 



5.5 Numerical Analysis 

Using the second statement of Theorem 5.3, in the presence of condensation we may obtain 
an expression for the short cycle density for a gauge-broken IRH Bose-Hubbard Model. 

Assuming our conjecture is correct, we apply some simple numerical techniques to this 
expression in order to compare long cycles with the Bose-Einstein condensate (see Appendix 
C for a brief explanation of the algorithms used and the numerical source code). As may 
be seen from Figures 5.3 and 5.4 the calculations certainly agree with Lemma (5.3), i.e. 
that the absence of condensation implies the lack of long cycles and visa versa. However 
more importantly they also indicate that while the presence of condensation coincides with 
the existence of long cycles, their respective densities are not necessarily equal. In fact, it is 
indicated that the long cycle density may be greater than or less than the condensate density 
for differing parameters. 
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Fi gure 5.4: Comparison of long cycles with ^^ onc j for various values of (3 and X 
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Penrose's Treatment of the 
Infinite- Range-Hopping 
Bose-Hubbard Model with 
Hard-Cores 



Here we shall briefly summarise the procedure in Penrose's work [66] in deriving an expression 
for the condensate density for the Infinite- Range- Hopping Bose-Hubbard model. We consider 
a collection of n bosons on a lattice of V sites with the hard-core condition, i.e. no more 
than one particle per site. As we are dealing with bosons, it is not necessary to specify which 
particle is where, only which sites are occupied and which are not. The allowed states of the 
system lie in a n-particle subspace H„ of tL^. Let X be the set of lattice points which are 
occupied. Then denote the corresponding quantum state \X) G T-L n . For a given n, there are 
( y ) such states, hence dim"H n = 

Define creation and annihilation operators A* and A by the following: 




and 
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These operators satisfy the commutation relation [A, A*] = 1 — yN, where N is the number 
operator defined by 

N\X) = #(X)\X). 

Penrose finds that the Hamiltonian for this model may be expressed in terms of these oper- 
ators 

H = N- A* A. 

where the A* A encompass the kinetic term of the Hamiltonian. 

From the commutation relation [N, A*] = A* one may interpret the operator A* as one 
which increases the number of particles in the system by one. However defining the following 
operator 

L = N — A* A - V~ l (N 2 -N) = H - V~ l (N 2 - N) (A.l) 

one may calculate that [L, A*] = 0. Note also that [H, L] = [H, N] — 0. Hence we have 
two simultaneously diagonalisable quantities with which we may be able to parameterise the 
eigenstates of H . By (A.l), we immediately have 

H = L + V~ l (N 2 - N). 

Let r be any integer satisfying < r < V/2 and construct the following subspace of "H: 

<H T = {\<P):N\<P)=r\<j>)} 1 

then by construction every state \<f>) 6 "H r satisfies 

N\<f>) = r|0), H\<f>) = r|0), L\<j>) = [r - V~\r 2 - r)]\<f>). 

Then consider the further sequence of subspaces 

7~Lr,m = (A*) m r T-i T 

for m — r, r + 1, . . . , V — r (with the consequence T-L r , r — H r )- One may then observe that 
for every state \<f>) G % r ,m safisfies 

N\<j>)=m\<t>) 1 L\<p) = [r-y _1 (r 2 -r)]|0) (A.2) 
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and consequently 

H\<f>) = [r- V~ l {r 2 -r) + V~\m 2 - m)]\</>). (A.3) 

The subspaces "H r . m for all < r < V/2 and m = r, r + 1, . . .V — r are all orthogonal, as 
they correspond to unique pairs of eigenvalues for the operators N and L. And thus they 
are eigenspaces of H. 

Penrose then determines the dimensions of these subspaces to be 

fV\ fV - 2r + 1 
dimH r , m = dimft r = I j I y_ r + 1 

with which we can write down an expression for the partition function of an n particle system 
as a scaled sum of contributions from each subspace % r ,n : 

min(n.V— n) 

Z p (n,V)= £ z(r,n 7 V,P), 

r=0 



where 



z{r,n,V,P)= ( ^y^+l ) (r) 6Xp {"f ^ Vr ~ r * + T + n * ~ ^ } " 



A.l Condensate Density 



To find the condensate density, we need to evaluate 

eS** = (A.4) 

V— >oo V 

where 

trace [A*Ae~ 0H ] 
^ N °' = trace [e^] 

Again we use the decomposition of the Hilbert space % n into the subspaces H r>m as performed 
above. By (A.l) we have that A* A = N — H, and using (A. 2) and (A.3) we see that "H r>m 
is an eigenspace of A* A with eigenvalue 

m — [r — V~ 1 (r 2 — r) + l/ _1 (m 2 — m)\ = — [(m — r)(V — m — r + 1)] . 
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This gives us sufficient information to fully diagonalise the traces in (A. 4) to obtain 

£cond = .l im 



/ -l im M 

icond — J" 11 \r 



mm(n,V— n) . x . . 

= 1 E (Af-rKV-^-r+l) 

To evaluate this expression, we compare with Penrose's Large Deviation-like variational 
formula from Section 4.4. Set x = r/V and p = N/V and obtain 

h v (x) = (p - x)(l - p - x + 1/V). 

This fa v converges uniformly in [0,min(p, 1 — p)] as V — > oo to a continuous function h : 
[0, min(p, 1 — p)] — y R of the form 

h(x) = (p — x)(l — p — x). 

Inserting this into (4.14), we obtain the density of the condensate for this model: 



^cond 



o if pe [o, P/3 ]u[i- P/3 ,i], 

{p-pp){l -p-pp) if pe[pp,l-pp] 



where pp is the critical temperature as stated in Chapter 4 (see equation (4.3) and the 
following paragraph for a definition). 
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Counter-example to the Positivity of 



In Chapter 3 a series of three conditions was developed which when satisfied by a particular 
model implies that the density of the Bose-Einstein condensate equals the long cycle density 
in that model, i.e. g£ = g^ ond - Here we wish to give an example of a model which breaks 
Condition I, which demanded that for any single-particle operator A > 0, the respective 
cycle expectation is non-negative, i.e. Cy(q,A) > 0. In this appendix we shall analyse more 
deeply the Infinite-Range-Hopping Bose-Hubbard Model with Hard-Cores of Chapter 4 in 
order to evaluate the cycle statistics of this model for the particular case of three particles 
hopping on a lattice of V sites. After explicitly finding the eigenvalues and eigenvectors for 
the three particle Hamiltonian if 3 c v (using the notation of Chapter 4), we derive expressions 
for projections upon this eigen-basis and rewrite the partition function with respect to this 
diagonalisation. Then we can obtain an analytical expression for the cycle density Cy(q) for 
any q. We find that the addition of a the hard-core interaction actually causes a damping 
effect on the cycle density of the particles. This damping can be exploited by choosing a 
particular operator which projects onto a vector perpendicular to the ground state, whose 
cycle expectation can be negative for certain temperatures. 
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B.l Diagonalising the 3-particle Hard-Core Hamiltonian 

To begin, recall that the one particle Hamiltonian h v has eigenvalues and 1. The eigenvector 
corresponding to is g = V^ _1 / 2 (l, 1, 1, ... 1) so then we may write h v = 1 — |g)(g|. Let 
ej = (0, . . . , 0, 1, 0, . . . , 0) with 1 in the 2-th place, which represents a particle on the site 
labelled % on the lattice. One may then evaluate that 

3 

Note also that for any normalised clg with components q, then 

5>c 3 - = <% - 1/K (B.l) 

C-Lg 

i?3 y is the Hamiltonian for three particles with Hilbert space T-V^y Any state <fi G n^ v may 
be expressed in the form 

v 

(f)= ^2 <fiijk &i ® ej <g> e fc (B.2) 

i,j,k=l 

where 0^ G C with the restriction that 0^ = when i,j,k are not distinct, i.e. <paj = 
4>iji — 4>ijj — 0, to satisfy the hard-core condition. Denote $ = J^ijk^ijk- 

We want to find such states (i.e. values of <f>ijk) that solve the eigen-equation H^ c v (p = Ecf>, 
where E are the energy eigenvalues of the system. Inserting (B.2) into the eigen-problem, 
one obtains 



Hly(p = 30 - — ^2 fcjk I ^2 e i> ) ® e j ® e fc - 7F e i ® ( e -?' ) ^ 

i,j,k=l \i'=£j,k / i,j,k=l \j'^:i,k / 

y Yl ^ k ^ ® ^ I S 



efc/ 

i,j,k=l \k'^i,j 



ej Cg> e 7 - ® e fc := 
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which implies that 



Hjk 



(3-E)V 



'ijk' 



i,j,k distinct, 
otherwise. 



(B.3) 



Suppose for the moment that E ^ 3, then for distinct indices k denote 



2 i fc ~~ Ys ^ v i k bki ~ ~Yj 



'k C-ij 



'ijk' 



(Note that we can drop the restrictions in the sums here since <pij k = when k are not 
distinct). Then by definition an = b u = Cu = 0, J2jk a jk = J2ki b ki = J2ij c ij = $ and 



J2 = Y bki Y bki = Y Ci i Y a > k = Y 



Cij. 



Using this notation, (B.3) becomes 



k — A 1 [ajk + b ki + c^) 



(B.4) 



where A = (3 - E)V . 
Sum (B.4) over j, k): 



22 'Pijk = A 1 22 ( a i k + bki + Ci i) 

i(¥=j,k) i&3,k) 



ajk — A 1 



{V-2)a jk + ^ b ki + ^2 c v 



(X — V + 2)a jk = Y bki + Y Ci i ~ bk 'J ~ ° k i- 



By symmetry, summing (B.4) over i, k) and over k(^ one can obtain the following 
system of equations 



(A - V + 2)a jk + b kj + c k j = /]bkj+)]ci 

i i 

a-ik + (A - V + 2)b ki + c ik = ^2 a Jk + Y c 

j j 

dji + bji + (A - V + 2)dj = ^2 a Jk + Y b 



(B.5) 
(B.6) 
(B.7) 
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We want to find all possible solutions to this set of equations, so that to obtain values for E 
and corresponding information on the coefficients 0^ (i.e. the c$j). Let 

Ak = Y^ o,jk B k = bjk Ck = Cjk 

j j j 

= Y a i k = Y b i k = Y c 'j k 

k k k 

where again J2k A k = J2k B k = Y,k C k = Aj = J2j Bj = J2j Q = Also we nave tliat 



Ak — B k , B k — Ck, Ck — A 



k- 



For brevity also denote — a = A — V + 2, then using these additional notations the simultan- 
eous equations (B.5)-(B.7) may be written as 



-aajk + bkj + Ckj — Ak + Cj 
&ik — &bki + Cik = Bk + Aj 



a ji + bji — &Cij = C k + Bj 



(B.8) 
(B.9) 
(B.10) 



which we may rewrite in matrix form as 





-a 











1 


1^ 




a jk 




(A k 


+ C> 







— a 





1 





1 




bjk 




B k 


+ Aj 










—a 


1 


1 







c jk 




c k 


+ Bj 







1 


1 


— a 










a kj 




Aj 


+ c k 




1 





1 





—o 







bkj 




Bj 


+ A k 


V 


1 


1 















{ckj 






+ B k/ 



(B.ll) 



Sum each row of (B.ll) over both indices 

Begin by summing the first row of (B.ll) (which is equation (B.8)) over both j and k (j ^ k), 
then 

-a ajk + Y hk i + Y Ck i = ( V " ^ Y bki + ( V ~ ^ Y Ci i 

U,k) (j,k) (j,k) i,k i,j 



-cr + 1 + 1)$ = 2(V - 1)$ 
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-a - 2(V - 2))$ = 



which implies either 



a = 2(2-V) or $ = 0. 
This result is independent of the row chosen, due to the symmetry of the system. 

Case I: a = 2(2 - V) =>- E = £ 

In this case we obtain the following for (B.ll) 





^2(V -2) 











1 


1 \ 








' A k + 







2(V-2) 





1 





1 




bjk 




B k + Aj 










2(V- 


2) 1 


1 







Cjk 




Ck + Bj 







1 


1 


2(V -2) 














Aj + C k 




1 





1 





2(V-2) 







bkj 




Bj + A k 




{ 1 


1 











2(V-2)j 








{Cj + Bkj 


which has solution 






















= bjk = cjk = 


O'kj = bkj = 


Ckj = 1 









since then Ai = Bi = d = V — 1. Therefore using (B.4), the corresponding normalised 
eigenstate for the energy eigenvalue 6/V is 



X, 



1 



V ^-i)(v^) J>^ 8e - 



This represents three bosons all in their lowest energy state. 
Case II: $ = where A^B^Ci are not all zero 

We shall sum each of the three simultaneous equations (B.8)-(B.10) over j). Taking 
the first equation we get 

-aA j + B J +Cj = S ^ j B k -Bj + (V -l)Cj Aj + Bj + aCj = 

k 

=<5>=0 



where a = —a — V + 2 = V(l — E) + 4. Performing the same to the other two equations 
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(B.9) and (B.10) one may form the system 



(s 


1 


A 


w 


1 


a 


1 






1 







0. 



(B.12) 



Since we are assuming Aj, Bj, Cj ^ then this system of equations is solvable if and only if 
the determinant of the matrix is zero, and therefore a must satisfy: 

(a- l) 2 (ff+2) = 
E = 1 + 3/V or £ = 1 + 6/V. 

Each of these eigenvalues will be considered separately. 

Case 11(a): $ = where A4, B { , Ci are not all zero — taking E = 1 + 3/V 

In this case, <r = 1, and hence (B.12) becomes 



i i i\ 

1 1 1 
1 1 ly 



(a \ 



A,- + B.- + C,- = 0. 



The matrix has rank 1 and thus has a two dimensional solution space, which implies that 
this eigenvalue is doubly-degenerate. Consider the following two orthogonal solutions 



Aj = -Bj with Cj = 0, 



and 



— Bj — \Cj 



and take them one at a time: 
• Aj = -Bj with Cj = 
For this possibility, a = 1 



a = 1 — V and so (B.ll) will become 



fv-i 











1 


1 \ 




djk 




f A ) 





V - 1 





1 





1 




bjk 




Aj — Ak 








V - 1 


1 


1 







Cjk 




-Aj 





1 


1 


V - 1 










a-kj 




Aj 


1 





1 





V -1 







b k j 




A k - Aj 


V 1 


1 











V-lj 




\Ckj) 




{ -A k ) 
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which has solution 



1 

V 



( A^{V-\)A h \ 
V-2 

Aj - A k 

V-2 
A k +{V-l)A 3 
V-2 

-(A, - A k ) 

A J+ (V-l)A k 
\ V-2 ) 



dj k C k j 
Cj k Cl k j 

bjk = —bkj 
and djk + bjk + Cj k = 



Therefore using (B.4), the corresponding eigenstate for the energy eigenvalue 1 + 3/V may 
be written as 



X 7« E 



Aj + (V- l)A k _ + Ak _ A% _ Aj + {V- 1)A 



V-2 



V-2 



e; ® e 7 - (g) e fc 



which when simplified and normalised (i.e. taking £\ |A;| 2 = 1 with '^2 i A i = 0) becomes 

1 



X 7 



_= (A k - Ai) e* ® Gj <g> e k . 



• Aj — Bj — \Cj 

By following the same procedure as above, one may find a second eigenvector of the degen- 
erate eigenvalue 1 + 3/V of the form 

1 



X, 



-=y E ( Ai + Ak ~ 2A ^ * ® ^ ® ^ 



which is orthogonal to X7. 

Case H(b): $ = where A», Bi, Ci are not all zero — taking E = 1 + 6/V 



In this case, a = —2, and hence 


(B.12) is 








(-2 1 








1 -2 


1 






I 1 1 


- 2 J 


w 



0. 
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The matrix has rank 2 and thus has a single dimensional solution space defined by: 



A? — Bj — Cj. 



As a = — 2 =^> cr = 4 — V and so (B.ll) will become 



V-4 











1 


1 


\ 






(l\ 





V-4 





1 





1 




b jk 




1 








V-4 


1 


1 







Cjk 




1 


















= (Aj + A k ) 







1 


1 


V-4 










O-kj 


1 


1 





1 





V-4 







hkj 




1 


V 1 


1 











V - 


V 


\CkjJ 




w 



which has solution 



Cljk bjk Cjk Qikj bfsj C-kj 

Hence we can construct the respective eigenvector of 1 + 6/V: 

1 



Aj + A k 



V 



This is a bosonic state of two particles in their ground states with one excited particle. 

Case III: $ = where A, = B { = d = 

In this case, (B.ll) reduces to the following eigensystem 



(° 











1 


A 



















1 





1 




bjk 




b jk 











1 


1 







Cjk 


= a\ 


Cjk 























1 


1 













a kj 




a kj 


1 





1 













b kj 




b kj 


V 1 


1 















{<*) 




\CkjJ 



where Iq is the 6x6 identity matrix. 
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The above matrix has the following eigenvalues and eigenvectors: 

a = 2 (1,1, 1,1,1, If 

a = -2 (-1,-1, -1,1,1, If 

a = 1 (1,-1, 0,-1,1, Of and (1, 1, -2, -1, -1, 2) T 

a = -1 (-1, 1, 0, -1, 1, Of and (1, 1, -2, 1, 1, -2) T . 

Note that the eigenvectors for the degenerate eigenvalues were chosen to be orthogonal. 
Considering each eigen-solution in turn: 

Case 111(a): a = 2 =>■ E = 2 + A/V 

Here (a jk , b jk , c jk , a kj , b kj , c kj ) T = (1, 1, 1, 1, 1, if so we have that 

ctj k bj k Cj k o, k j b k j c k j 
and hence inserting these components into (B.4) we can construct a normalised state 
X 5 = = (a jk + a ki + a,ij) e { (g e,- <g> e k 

where = dji and Aj = implies that a^- = for all j — 1 . . . V, \ a ij\ 2 = 1- This 
is a bosonic state of two excited particles and one in its ground state. 

Case 111(b): a = -2 E = 2 

Here by (— 1, — 1, — 1, 1, 1, if one deduces the relations aj k = —a k j = bj k = —b k j = Cj k = 
—c k j with which the normalised state can be written as 

Xi = —= ^ { a jk + a-ki + dij) Gi <g e,- <g e k (B.13) 

V o V ... 

where now = —aji with Yli a ij = f° r an J = 1 • • • V", Ylij \ a ij\ 2 = 1- This may be 
rewritten in a more intuitive form by setting to obtain 

Xi = (g <g Ci <g c 2 - g <g> C 2 g) Ci + C 2 <g> g <g Ci - Ci <g g (g) C 2 + Ci (g) c 2 <g g - C 2 <g> Ci <g g) 

V6 

where Ci _L c 2 _L g. This represents a fermionic state with a single particle in the ground 
state. 
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Case III(c): a = 1 => E = 2 + 3/V 

This eigenvalue is doubly degenerate, so we shall consider each eigenvector separately: 

• By (1, —1, 0, —1, 1, 0) T , have that dj k = —a k j = —bj k = b k j and Cj k = c k j = and hence 

Xn = 1 ^ (ojfc - a ki ) e { <g> e 3 - g> e k 
V 2 ( V ~ 3 ) {i:j , k) 

with a,ij = —dji where the scaled as above. 

• By (1, 1, -2, -1,-1, 2) T , a jk = -a kj = -b jk = b kj = \c jk = -\c kj and hence 

X12 = - 1 V] (a jfc + a fei - 2^-) e* <g> ® e fc . 
V " 3 ) (S) 

Case 111(d): a = -1 =>- £7 = 2 + l/V 

• By (—1, 1, 0, —1, 1, 0) T , aj k = a k j = —bj k = —b k j and Cj k = c k j = and thus 

X 9 = 1 ^ - Ofti) ej ® e^ <g> e k 
V 2 y V ~ l ) (i,j, k ) 

with and again the a^-'s scaled as above. 

• By (1, 1, -2, 1, 1, -2) T , a jk = a kj = b jk = b kj = -\c jk = -\c kj to form 

X] + afei ~ 2ai ^ e * ® e J ® e?c - 

Case IV: £ = 3 

Here we cannot apply the above method as the right- hand- side of (B.4) is not defined. 
However we can propose a fermionic state with all excited particles of the form: 



X 2 = —= y~] sgn(vr)c 7r( i ) ® c n{2 ) <8> 
V 6 



(3) 



where Ci _L c 2 _L c 3 _L g and ||cj|| = 1. This can be quickly verified by recalling that 
h v — 1 — |g)(g| and H^° v = V^ c hyV^ c . Similarly for bosons: 

X 6 OC ^l<3<3 J] ® ^-(2) ® ^-(3) 

where ^ = 0^- = 0^- = 0^ = and ^ = J2j 4>ijk = J2k fart = °- 
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APPENDIX B. COUNTER-EXAMPLE TO THE POSITIVITY OF C(Q,A) 
Summary of Eigenvalues and Eigenstates of H^ c v 

We split the Hilbert space into a symmetric subspace, an anti-symmetric subspace, and 
a subspace which is neither symmetric nor anti-symmetric. On the odd (fermionic) space 
the three-particle Hamiltonian H% c v has eigenvalues 2 and 3 with eigenvectors: 
• eigenvalue 2: 

Xi = — = ( a jk + a ki + e, ® ej <g) e k 



where a^- = — aji with = 0, £\ • |ajj| 2 = 1 and ^ c%- = for all j. • eigenvalue 3: 

1 

where Ci _L c 2 _L c 3 _L g and ||cj|| = 1. 

The eigenvectors of c v in the even (bosonic) space are of the form: 

• eigenvalue 6/V: 

= — , e, (g) e,- cg> ei,. 

^V(V-l)(V-2)^- 1 

• eigenvalue 1 + 6/V: 

X 4 = ■===, y~] (Ai + Aj + A k ) e { ® ® e fc 

V3(F-2)(T/-3) ( ^ V 

where £.4, = and EJA| 2 = 1- 

• eigenvalue 2 + 4/V: 

X 5 = 1 (a jk + a ki + ay) &i ® e,- <g> e fc 

where a^- = with an = 0, \ a ij\ 2 = 1 an d J2i a ij = f° r an J- 

• eigenvalue 3: 



X 6 0C ^2*3 Yl e V(D ® ^.(2) ® e M3) 

<*li«2,«3> 7T6S3 

where Cjjj ^jjj ^ and Cjjfc ^ ' • Cj^fc ^ijfc 0- 
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The following eigenvalues of H% c v are degenerate with eigenvectors that lie outside both the 
even and odd eigenspaces: 

• eigenvalue 1 + 3/V: 

X 8 = - 1 V (A4 + A k - 2Aj) ei <g> ej ® e k 

VQV(V - 2) ^ 

where = and £\ \ A i\ 2 = 1 - 

• eigenvalue 2 + l/V: 
X 9 = 7 = ^ ( a jk ~ a ki) £i ® (8) e fe 



where = a^ with an = 0, Ylij \ a ij\ 2 = 1 an d Sj a «j = f° r au i- 



= ^ (ojfe + a fci - 2ay) e» <g> ej (8> e fc 



• eigenvalue 2 + 3/V: 

Xn = /o / T ! ^ Yl ( a i k ~ a ^ ei ® e J ® efe 

X12 = , ^ (oj-jt + a ifc - 2ay) g) ej (8) e fc 

V 6 ( V - 3 ) <i,j,fc> 

where = —aji with an = 0, ^ . | 1 2 = 1 and X^ a «j = f° r an i- 

Note that these eigenvalues can also be acquired from equation (A. 3) by Penrose (see Ap- 
pendix A): 

r — V~ l {r 2 — r) + V~ l {m 2 — m) 
where here m = 3 and r = 0, 1, 2, 3. 



We shall next derive explicit expressions for projections upon these eigenvectors. However 
in order to do this, one must consider the two particle case for some expressions which will 
be required for the three particle calculation. 
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B.2 The 2-particle Hard-Core Hamiltonian 

By performing the technique of the previous section for a system of two hard-core particles 
on a lattice of V sites, one can diagonalise the Hamiltonian H^y as follows: 

In the odd space the two-particle Hamiltonian has eigenvalues 1 and 2 with eigenvectors: 
• eigenvalue 1: 

1 V 1 

Yi = -j= (ci - cj) ei <8> = — = (g g c - c g g) 



where c J_ g and ||c|| = 1. 
• eigenvalue 2: 



1 



Y 2 = - ciij ei g ej = "7=(c g c' - c' £ 
2 . y2 



where ajj = — with = 0, £\ . \dij\ 2 = 1 and ^ = for all j. The second expression 
is found by taking a^j = 2 _1 / 2 (qCj- — c^Cj), in which case c, c' J_ g and clc' with unit norm 



The eigenvectors in the even space are of the form: 
• eigenvalue 2/V: 



1 V 1 / V 

Y 3= .^ T _ 1) =E( 1 -^) e ^ e ^ = ^/^=Y ) ^g®g-E e ^ e " 



eigenvalue 1 + 2/V: 



v 

Y 4 



^(V-2) ^ ^ + ~ ^ ^ ® ^ 



1 ( v 7 ! 7 (g ® c + c g g) - 2 a e; g e 4 j 



where c _L g, ||c|| = 1. Note that X 4 (c) J_ X 4 (c') if c J_ c'. 
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• eigenvalue 2: 

I v 

Y 5 = — = 22 a iMi ® e 'j + e i ® e i) 

where a i3 - = with = 0, Yn a ij = V j and Y%j=i \ a v\ 2 = 1- 

Denote P/ as the projection onto Yj. We may evaluate these five projections as follows: 

c:£ci=0 c _L g 



1 

2V 



using (B.l). 



4 

a c,c' 

C-Lc'-Lg 

-[-D(ii,«2,ii,j2) + D(i 2 ,ii,j 2 ,ji) ~ D(k,i2,32, ji) ~ D(i 2 ,ix,j 1 ,j 2 

(B.14) 



where 



j 2 Ci 2 Cj 2 



D(i 1} i 2t 31,32) = Y 44= Y Ci ^[ £ 4 c j 

c: c_Lg c':c'_Lc,g c:c_Lg \c':c'_Lg 

^ ' ^ ' c:c_Lg 

When we reinsert this into (B.14), the latter terms of the above expression will all cancel 
leaving 

p\ (z'i, i 2 ; ji, j 2 ) = ~(Si 1 j 1 5i 2 j 2 — Si 1 j 2 6i 2 j 1 ) — —(5^! — Si 1 j 2 — 8^ + 5i 2 j 2 )- (B.15) 
Similarly we have 



^4 (*1>*2; J'1^2) = 77777 ^(<W + + <W + $aja _ 4 /^ / )( 1 ~ *tii a )(l ~ ^'ij 2 )- 



2(y-2) 

We cannot derive this projection from the eigenvector Y5 itself due to lack of information, but 
since Pf ] = I® - Pf ] - P 2 [2] - Pf ] - Pf ] and the matrix for I® has kernel 8 hh 8 i2h {l - 6 t 
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we can deduce 



+ (v _ 1) 1 (v _ 2) |(l-^ 2 )(l-^ 2 )- (B.16) 

Consequently the partition function for the two particle hard-core Hamiltonian may be 
written as 

e -^v = e -/»pM + e-^Pi 21 + e^Pj 21 + e -^/v)pM + e -VpW 

= e-20/Pl + ( e -/3 _ e -^)pf ] + (e" 2 ^ - e- 2/3 )Pj 21 + (e"^ 1+2 /^ - ^)pf ] . 

With this one can evaluate, for instance, the 2-cycle density 

v v 
c 2 (g = 2) = c *(hi) where c ^h3) = ^2e- pH ^v(i,k;k,j) = f,(V) + g,(V)6 ir 

i,j=l k=l 

Observe that the term g /3 (V)8ij is a consequence of the hard-core interaction, so deleting it 
returns the cycle density of the two free particles. It can be shown that for two particles the 
functions f and g p are both non-negative for all possible cycle lengths (i.e. one and two). 
However in the three particle case, the Hilbert space can be divided into subspaces which 
are not only even and odd, but also those which are neither even nor odd. These additional 
spaces can result in negative contributions to the cycle density expression, as we shall shortly 
see. 

B.3 Cycle Density Expression for the 3-particle Hard-Core Model 

Returning to the three particle case, let Pj denote the projection onto the eigenvector Xj. It 
is useful to set the operator 

A ± (i, k; j, I) = Sij + 5 k i ± (S u + 5 jk ) 
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and also to abuse notation slightly by denning 



Si i 



ijk 



1 iii^j^k 

otherwise. 



Taking these projection in order of complexity. First one may easily see that 



1 

v(v - rj(v - 2) 



Pa(ii,i2,ia'Ji,j2,j3) = ttttj — ttttt — ^(1 - S ili2la )(l - S jlhh ) 



where we introduce 1 — terms to enforce the summation's conditions (and hence the 

hard-core condition) that H- 

Secondly, the eigenvectors X4, X7 and X§ all contain terms of the form Ai where Yli Ai = 
and Yli \Ai\ 2 = 1. Then choosing X4 we have 

P±{iui2ih\3\i32,h) = X] ( A h+ A i 2 + A i 3 )(Aj 1 + A h + A h ). 

But since A4 = one can say that the vector A _L g (as we do in the two-particle case) 
and hence 

^ Ai-Aj = Sij - —. 

A±g 

With this we may calculate the following (with Pi m := Pi + P m ): 

9 



P4{il,i2,i3;jl,32,j3j 



3(V-2)(V-3) 



V 



^7,8(«1^2^3; 31,32,33) 



x (1 &iii2h)(X fijij2js) 

1 



A (ii,i 2 ; 31,32) + A fa, is; 32,3s) + A (i 3 , h; j 3 , j x ) 



3V(V -2 

X (1 — ^1*213) (1 ~~ ^jrj2h)- 



Thirdly consider the eigenvectors Xi, Xu and X12. These all have the property that = 
— dji. From the expression (B.13) for the eigenvector Xi we may construct the corresponding 
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projection P\: 

Pi(h, k, 32,33) = gTF '^2(a i2 i a + a iah + a ili2 )(a j2h + a jajl + a jlh ). 

a 

In two particle calculation above for odd we derived a projection (B.15) of the form: 

P\ (i\,i 2l jl) 32) = ^ ] a i\i-fihh = 2 (^Ul^aj^ ~~ ^«lj2^2jl) ~~ l^yi^ixjl ~ ^hj2 ~ ^231 + ^2^2) 
a 

which we can use here to obtain 



Pi(ii,i2,h]ji,j2,j3) = 7^7 



rol [9] [ol 

P-2 (i2,k;32,j3) +P2 («2,«3;i3,jl) + ^2 (*2, is] jl, 32) 

[21 \2\ \2\ 

+ P 2 (i3,ii;j2,33j + P 2 («3,n;i3,ii) + ^2 {k,k; 31,32) 

[21 T2l [21 

+ P\ («l,«2;j2,j3) + ^2 (il,i2; 33,jl) + *2 {k,i2,3l,32) 

This can be simplified (all the 7^7 terms cancel) to obtain 
Pi{ii,k, 13,31,32,33) = ^yhhjA~(i 2 ,i3; j2,j3)+^ 1 i2 A ~(^3,^; ji, j3)+S ilj3 A~(i 3 ,i 2 ; j 2 ,ji) 

+ fii2jl(,fi'hj2 ~ ^3i3) ^«2i2 (^«3J3 _ ^«3il) ^2.73 (^3.71 _ ^«3J2 < 

In a similar fashion, one can find expressions for the projections Pu(ii, %2, hlji, 32, 33) and 
Pviiii, h, h;3i, 32, 33), which we shall combine to obtain 

Pu,i2{ii,i2,i3,3i,32,33) = Q(y 1 _ 3) j^m A + (i 2 ,i 3 ; j 2 , j 3 ) + 5 i2h A + (ix,i 3 ; j x ,j 3 ) 

+<5j 3J3 A + (z 1 ,z2; ji, 32) + ^ 3 j2( — 25 i2 j 3 — 5 ii;7 - 3 — 5j 2jl ) + ^gj^— 25j i: ,- 3 — 5j 1J2 — 5 i2 j 3 ) 

3 

-25i 2jl 5i lj2 -5i 2j3 5 ilj2 -S i2jl 6 ilh - —[A~(i 1 ,i 2 ;j 1 ,j 2 ) + A~(i 2 ,i 3 ;j 2 ,j 3 ) 



+ A 0' 3 ,n; j 3 , jx)\ |(1 - <Ji 1 i ai3 )(l - S hj2j3 ). 

Finally the remaining vectors we shall consider, X5, X9 and X10, all have that = a^. As 
above, these three-particle projections can be expressed in terms of a sum of a particular two- 
particle projection with different arguments. In this case the projection is P 2 %2\j\, J2) 
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for which we have an explicit expression (B.16). With this one can evaluate 

1 f 9 2 

Pbih^MjuhJa) = s(V -4) \ {V -1){V -2) ~ V^2 ^ hjl + 6ilh + " ' + Siaj ^ 



1 

+ 2 



5 iljl A + (i 2 ,i 3 ;j 2 ,j 3 ) + S il j 2 A + (i2,i 3 ;j 1 ,j3) + 6 ilj3 A + (i 2 ,i 3 ;j 1 ,j2) 



J f^i2ji\^i3j2 ^iajs) ^252 (^353 ^351) ^253(^351 ^352 1 



,1 ^112*3) (^ ^515253) 



^9,10(^1,^,^; 31,32,33) = Q(y _ j^l51 A («2,«3; J2,J3) + ^252 A ^ Jl, J3) 

H-^jgJgA («l,22j Jl? + ^352 (2^253 _ ^153 _ ^ji) ^351(^^153 _ ^152 _ ^253) 

+28 i2jl 5 ilj2 -5i 2j3 5 ilj2 - 5 i2jl 8 ilj3 - y _ 2 [A~(ii,i 2 ; ji,j 2 ) + A~(z 2 , z 3 ; J2, J3) 



+ A (z 3 ,Zi; J 3 ,Ji)] |(1 - <Ji li2ia )(l - (Sjijaja). 

What remains is P 2j 6, the projection onto the space spanned by X 2 and Xg. This can be 
written as P 2 ,6 = I ~ Pi — P3 — Pa — P5 — Pi,& — -^9,10 ~ Pn,i2 where we note that / can be 
written as 5i 1 j l 5i 2 j 2 5i 3 j 3 (l — 5^^). 

The Hamiltonian H$y is diagonal with respect to these vectors X, and hence it can be 
written in terms of these induced projections: 

e -^,V = e-^Pr + e~^P 2fi + e~^ v P 3 + e~^ 1+6 ^P, + e" 2/3 P 5 + e^+^P^ 
+ e -/3( 2+ i/v)p 9il0 + e -^ + 3/y) Piiii2 

= a" 3 / 3 / + (e- 2 ^ - e-'^P, + (e" 6/3 / y - e" 3/3 )P 3 + (e^ 1+6 ^ - e" 3/3 )P 4 
+ (e-2/3 _ e -3/3)p 5 + (e -/3(i+3/y) _ e -^)P 7j8 + ^-P^m - e - 3/3 )P 9 , 10 

+ (e"^ 2+3 / y ) - e- 3/3 )Pn,i2. 

Thus we can find an explicit expression for the kernel of the operator e~^ H 3,v _ L e ^ us consider 
the case of a cycle of length three, i.e. 

(k,i) 
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Then one can evaluate the three-cycled projections to write 

£)j(z,M;M,j) = 

(k,i) 

E P l(j , M;M , i) = o^^) + ^ 



(k,i) 



J2P 3 (i,k,l;k,l,j) = 3 + <E - — 



U(y-l) 1) 



(k,i) 



VRfi fe f-fe I ;) ~ 3 )(^ ~ 7 ) ■ j (y-3)(y + 2) 



(k,i) 



C 1,1 1,1 ^ 2(V-3) V + 3 



(k,i) 



^ 9,10V, , , , 6(y-i) J 6(y-i) 



(fc,i> 



2_^Pii,i2{i,k,l;k,l,j) = — dij— 



6V J 6V 



(k,i) 

which returns the following expression 

^e-^3 hc v(,, k} i- k , I J) = f p (V) + 5 ij9f) (V) 
(k,i) 

where fp{V) is a non-negative expression of the form 
UV) = (e-^-e^)^^^+(e-^ 



-2/3 _ r -3P\ ( V 3 )(^ 7 ) _ f p -/3(l+3/V) _ p -3/3\ 2 ( V 3 ) 

j 6(y-i) 1 J 3U 

-^(2+i/v) _ r -^ (V-3)(V-4) _ , -^ (2+ 3/y) _ r - 3 ^ (^-2)(^-6) 
J 6(V-1) 1 J 6U 
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and Qpiy) is a not-always-positive function as may be seen from Figure B.l, 

V-2 , 



9p{V) 



+ [e -W/V _ e -3/S> 



v(v-r 



_|_ ( e -/3(l+6/V) _ e -3/3> 



3V 



+ ( e -2/3 _ e -3/3 



(V-3)(y + 2) _ ( , {1+3/v) _ 3P) V + S 
] 3(V-1) 1 j 3U 



j- e -/3(2+l/K) _ e -3/3 



(2^+l)(y-3) 
6(V- 1) 



,-/}(2+3/V) _ e -3/3- 



(2U + 3)(V-2) 
6V 7 



0.0001 H 




-0.0008 



Figure B.l: Plot of the function gp(V) versus (3 for V = 100 

Recall from equation (2.15) that Cy(q,A) is an expression for the expectation of a single 
particle operator A > in terms of cycle statistics. In the hard-core model, the expansion 
for three particles on a lattice of V sites can be written as (see (4.18)) 

1 



Z,(3,V)V 



trace 



7 , 3 C («A ) ®^A+ 9) ) 



(A®I^)U q B-^v 
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for q — 1, 2, 3. Taking the case q = 3, in terms of the kernel of H% c v , we may expand this to 
obtain 



4(3, A) 



Z,(3,V)V tTaCe ^ 



(A <g> I ® J) t/ 3 e-^ 



Z,(3,V)V 



(A (8) / ® /) t^e"^ 



(j,k,l;j,k,l) 



^A(j,m)e- pH ^(m,k t l-,k,l t j) 



Z,(3,V)V 



(j,k,l) m=l 

r 



^ ' ' 7,m=l 



Let A be a projection onto states orthogonal to the ground state g, for example let c = 
2" 1 /2( 1) _ l7 o,...,0) and fix A = |c)(c|. It has kernel A(j,m) = 2~ 1 (5 lj - S 2j )(S lm - S 2m ). 
Then 

v 



4(3,|c)(c|) 



1 



22,(3, V)V 



j,m=l 



g (v) 

Z P (3,V)V 



which is not necessarily a non-negative number (see Figure B.l). 
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Numerical Source Code 



Contained within this appendix is the C++ source code to generate Figures 4.1, 4.6 and 
5.4 of Chapter 5 on the Infinite-Range-Hopping Bose-Hubbard Model, which compare the 
density of the Bose condensate with the long cycle density. The algorithm represents the 
Approximated Hamiltonian, equation (5.2), as a finite dimensional symmetric tridiagonal 
matrix with respect to the occupation number basis, and diagonalises the matrix to find its 
eigenvalues and eigenvectors. 

The matrix is diagonal except in the case of r > 0, i.e. when BEC occurs. So for each 
H, one needs to calculate the solution of p'(r) = 2r, done using a simple binary search- 
like convergence algorithm (see lines 143-210). This value of r is then used to evaluate the 
pressure, density and the cycle-densities. 

This code uses the open source GMP [46] and CLN [48] libraries to boost the accuracy of the 
calculations beyond that of standard double precision arithmetic. This was done for two 
reasons, first so that the algorithm was capable of summing the exponents of the eigenvalues 
without overflowing, and secondly to zoom in on the critical values of fi to check for the 
discontinuity in r. The matrix diagonalisation algorithm is a modified tqli method from 
Numerical Recipes in C [69], enhanced for higher precision. The code automatically generates 
a plotdata file and a command file and calls gnuplot [96] to render two plots automatically, 
one comparing short cycles with the total density, the other long cycles and BEC. Both these 
plots are saved in two appropriately named EPS files. 
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i /* 

* soft-core-cycles.cpp 
* 

* To test hypothesis that \ sum_ { q=l} ~ \ inf ty c_V~\mu(q) = \rho (\beta , \mu) - r_0~2 

* And so r_0~2 is the density of BEC in the model 
6 */ 



#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 



<iostream> 
<f stream > 
<sstream> 
< string > 

< except ion > 

< iomanip > 
<vector > 
<cstdlib > 

< cmath > 
<cln / cln . h> 



using namespace std ; 
using namespace cln; 

/ * To enable Debugging output, activate the line below */ 
//Xdefine DEBUG 



/* Common functions */ 
26 cl_F f(int, double); 

cl_F f _approx ( int , int , double); 

cl_F part it i on_of _appr ox ( double , double); 

cl_F cycle_density_numerator ( int , double, double); 

double density (double , double); 
31 void der iv_of _pre s sure ( double , double, cl_F*, cl_F*) ; 

void tql i ( vector < cl_F >& , vector <cl_F >& , int, vector <cl_F >& , bool); 

cl_F pythag ( const cl_F , const cl_F); 



/* Advanced variables for precision purposes */ 



const int decimal_places = 40; 
const cl_F TINY_T = " 1 . e -40 _30 " : 
const double TINY = 1.0e-15; 



const int datapoints 



100; 



//How many decimal places the arithmetic uses 
//(as algorithm is arbitrary precision) 
//What constitutes negligibly small 
//Check " deciman_pl aces " before setting 
//set how accurate r value should be, e-12 is 
//good enough for 7-8 decimal places 
//Number of datapoints to take between mu_min 
//and mu_ max 



/* Initial defaults , usually overridden via command line */ 
double beta = 100; 
double mu_min = 0; 
double mu_max = 10; 
double lambda = 3; 

/* No need to change anything from here on */ 

f 1 o at _f ormat _ t precision = f 1 oat _f ormat ( de c imal_place s ) ; 

const cl_F zero = cl_f loat (0 . , precision); 



56 /* Handy macros */ 

#define SIGN(a.b) ((b)<zero ? -abs(a) : abs(a)) 
#ifdef DEBUG 

#define verbose if (0) ; else cerr 
#else 

61 #define verbose if (1) ; else cerr 
#endif 



//The matrix will be a symmetric tridiagonal . This is most efficiently stored in two arrays, 
//one with the diagonal entries , one with the of f - diagonal entries . 
66 vector <cl_F> diag; //tqli will place eigenvalues in this later 

vector <cl_F> offdiag; 
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vector <cl_F> eigenvectors; 
int number.of .occupations = 



//to store generated eigenvalues of matrix 



int 
{ 



main (int argc , char * const argv [] ) 

cl_F cycles=zero, expectation =zero , partition =zero ; 

double r _f ound =0 . , diff=0.0, dens i t y _f ound=0 . , shor t _ cycle s =0 . ; 

char barspin[4] = { ' \\ ' , 'I', '/', '-'}; 

int iteration=0, bariterator =0 ; 

double mu_iterator =0 . 1 ; 

of stream plotdata, gnuplot_commands ; 

//Check command line variables before starting 
if ( argc != 6 ){ 

cout << endl << "Usage: " << argv [0] 

<< " [number of occupations] [beta] [lambda] [mu_min] [mu_max]" << endl; 
return 1 ; 

} 

number.of .occupations = atoi ( argv [1] ) ; 
if ( number.of .occupations < ){ 

cout << endl << "Usage: " << argv [0] 

<< " [number of occupations] [beta] [lambda] [mu_min] [mu_max]" << endl 
<< "where you must supply a positive integer for occupations" << endl; 
return 1; 

} 



beta = atof ( argv [2] ) ; 
if ( beta < 0){ 

cout << endl << "Usage: " << argv [0] 

<< " [number of occupations] [beta] [lambda] [mu_min] [mu_max]' 
<< "where you must supply a positive integer for beta" << endl; 
return 1 ; 



<< endl 



lambda = atof ( argv [3] ) ; 
if ( lambda < 0){ 

cout << endl << "Usage: " << argv [0] 

<< " [number of occupations] [beta] [lambda] [mu_min] [mu_max]" << endl 
<< "where you must supply a positive number for lambda" << endl; 
return 1; 



mu_min = atof ( argv [4] ) ; 
mu_max = atof ( argv [5] ) ; 
if ( mu_min > mu_max ) { 

cout << endl << "Usage: " << argv [0] 

<< " [number of occupations] [beta] [lambda] [mu_min] [mu_max]" << endl 
<< "where mu_min <= mu_max ! " << endl; 
return 1; 

} 

//Construct filename for plotdata 
stringstream s; 

s << " plot data - bet a " << beta << "-lambda" << lambda << "-mu" 

<< mu_min << '"' << mu_max << "-n" << number.of .occupations ; 
string filename = s.strO; 

//Open "plotdata" file and prep are it 
plotdata.open(filename . c.str () , ios : : out) ; 

plotdata << "# Dataset for beta=" << beta << ", lambda=" << lambda 

<< " with " << number _ of _ oc cupat i ons << " occupations" << endl 
<< "# Data format: mu density r shortcycles " << endl; 



//Initialize the matrix 
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136 di ag . resize (number.of .occupations) ; 

offdiag . resize (number_of .occupations ) ; 

eigenvectors. resize ( number _ of _ occupations*number_of .occupations) ; 

//Set mu_iterator so that we evaluate correct number of datapoints 
141 mu.iterator = ( mu_max -mu_min ) / dat apo int s ; 

for(double mu=mu_min ; mu <= mu_max ; mu+=mu_iterator ) { 
bariterator=0; 
cl_F partition = zero; 
146 cl_F expectation = zero; 

#ifndef DEBUG 

//Progress meter 

iteration = f loor ( ( ( mu -mu_min ) * 100 ) / ( mu_max -mu_min ) ) ; 
cout << 1 \r 1 << iteration << "'/, complete"; 
151 cout . flush () ; 

#endif 

//Check to see if p ' (r) will ever intersect 2*r . Do this by evaluating both at TINY 
//and if p'(TINY) > 2*TINY => r>0 is a solution. Otherwise no. This is due to 
//concavity of p ' (r) as postulated (not proven, but looks likely). 

156 

r_found = 0.0; //set to lowest guaranteed solution 
//Find expectation of (a + a*) i.e. p ' (r) 

der iv_ of _pre s sure (0 . 0000 1 , mu , ^partition, ^expectation ) ; 
161 if ( expectation > 0.00002 ){ 

//otherwise skip mu loop to next iteration , since r=0 is only possibility 

//Now know r>0. Need to iterate over r to locate it. 
double r_increment = 0.2; //Set initial r_increment for r loop 
166 double r_initial = 0.002; //Set initial point of interval to seek for r in 

double r_final = 20; //Set final point of interval to seek r in 

// For this mu , find the max value of r s.t. 2*r = p ' (r) 

for (double r _ val =r _ ini t i al ; r_val < r_final; r_val +=r_increment ) { 
171 #ifndef DEBUG 

cout << ' \r ' << iteration << "'/, complete " << bar spin [ ( bar it er at or ++) °/.4] ; 
cout . f lush ( ) ; 

#endif 

//Find expectation of (a+a*) i.e. p ' (r) 

176 deriv_of .pressure (r_val , mu , &partition, ^expectation ) ; 

diff = double_appr ox ( expect at i on ) - 2*r_val; //difference 

if ( f abs ( dif f ) <= TINY ){ //if exp ect at i on==2*r , set r_found and break loop 
181 r_found = r.val ; 

verbose << ' ! ' ; 
verbose . flush () ; 
break ; 

} 

186 else{ 

//Has p ' (r) crossed 2*r yet? 

if ( diff > ){ //No, so leave all else as is 

if ( r.increient < 1.0e-15){ //Unfortunate possibility 
r_found = r.val ; 

191 break; //Almost gone beyond precision 

} 

verbose << 1 + ' ; 
verbose . flush () ; 

} 

196 else{ //Yes we've intersected the line 2*r , somewhere in interval 

// (r_val -r_int erval , r_v a I + r_ int erva I) 
r_final = r_ val + r_ increment ; 

r_val -= r_ increment; //need to double back , iterate & converge 

//to value of r at intersection 

201 r_initial = r_val ; 

r_increment /= 2; //use binary search as most reliable 

verbose << "->" << mu << " diff: 11 << diff << ", r_val : 11 << r val 
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<< " , r_increment : " << r_increment << endl ; 
verbose . flush () ; 



y//r loop 
J // end if (i . e . no intersection) 

verbose << " Found for mu=" << mu << " \ t r=" << r_f ound << endl; 

//Have found value for supremum of r (r_found), evaluate the pressure at r 

cycles = zero ; 

try{ 

//Eval sum of cycles - short cycle sum 
partition = partition_of _approx (mu , r_f ound) ; 

for(int q=l ; q<=number_of .occupations ; q + + ) { 

verbose << q << 1 \t ' << cycle_density_numerator (q , mu , r_f ound) 

/ partition << endl; 
cycles = cycles + cycle_density_numerator(q, mu, r_f ound ) ; 

> 

cycles = cycles / partition; // cycles /partition = cy cl e_ dens i ty 
//print values to file 

density_f ound = densit y ( r _f ound , mu) ; 
short_cycles = double_approx ( cycles ); 

plotdata << mu << ' \t 1 << dens ity_f ound << ' \t 1 << r_found << 1 \t ' << 
short_cycles << endl; 
} catch ( int a ){ //Divide by zero error, report, skip this iteration & continue 

verbose << "Divide by zero in TQLI , skipping evaluation for mu=" << mu << endl; 
continue ; 

} 

}// end mu lo op 

cout << 1 \r ' << "100*/, complete" << endl; //just tidy the output 
plotdata . close () ; 

// Generate a p I o t _ command file for gnuplot so graphs will be generated on the data found 
s . str ( " " ) ; 

s << " plot _ commands - bet a " << beta << "-lambda" << lambda << "-mu" 

<< mu_min << << mu_max << "-n" << number_of .occupations ; 

string plot_cmd = s.strO; 

gnuplot _ commands . open(plot_cmd . c_str () , ios : : out) ; 

gnuplot_commands << "set term postscript eps enhanced colour" << endl; 
s . str ("") ; 

s << " dens ity -r " 2 _ vs_short _cy cle s -bet a " << beta << "-lambda" << lambda << "-mu" 
<< mu_min << ' ~ ' << mu_max << "-n" << nuiber.of .occupations ; 

//Print plot_command file 

gnuplot_commands << "set out \"" << s.strO << ".eps\"\n" 
"set key spacing 1.3\n" 
" set key box\n " 

"set style line 5 It 1 lw 2\n" 
"set style line 6 It 3 lw 2\n" 

"set title \"{/Symbol b}=" << beta << ", {/Symbol 1}=" << lambda 
<< " with " << number.of .occupations << " occupations \"\n" 
"set xrange [" << mu.min << ":" << mu_max << "]\n" 
"set xlabel \"{/Symbol m}\"\n" 
"set key right bottom\n" 

"plot V " << filename << "\" using 1 : ( $2 - ( $3 * $3 ) ) title '{/Symbol r}-r*2' " 
"with lines linestyle 5, \\\n\"" << filename << "\" using 1:4 title " 
" '{/Symbol r}_{short}' with lines linestyle 6" << endl; 

s . str ("") ; //reset s for a new filename 

s << " r ~2 _ vs_long_cy cles -bet a " << beta << "-lambda" << lambda << "-mu" 
<< mu.min << '"' << mu.max << "-n" << number.of .occupations ; 
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271 gnuplot_commands << "set out \" " << s . str () << ".eps\"\n" 

"set key left top\n" 

"set title \"{/Symbol b}=" << beta << ", {/Symbol 1}=" << lambda 
<< " with " << number_of .occupations << " occupations \"\n" 

"plot \"" << filename << "\" using l:($2-$4) title '{/Symbol r}_{long}' " 
276 "with lines linestyle 5, \\\n\"" << filename << "\" using 1:($3*$3) title " 

IM r"2 1 with lines linestyle 6;" << endl ; 

gnuplot_commands . close () ; 
plot_cmd = "gnuplot " + plot_cmd; 
281 sy st em ( plot _cmd . c_ str ( ) ) ; //run gnuplot and point it to the p I o t _ command file. 

return 0; 

} 



286 

/* 

* Evaluates d/dr p(r) 

*/ 

void der iv_of _pre s sure ( double r, double mu , cl_F *partition, cl_F * expect at i on ) 

291 { 

cl_F temp = zero; 
♦partition = zero; 

for(int n=0; n<number_of .occupations ; n++) { 
296 diag[n] = f(n, mu) ; 

if (n ! =0) //note: offdiag[0] not used 

offdiag [n] = cl_float (r*beta,precision)*sqrt(cl_float (n, precision)) ; 
//Create identity matrix here, needed for tqli to work correctly 
for(int m=0; i<number_of .occupations ; m++) { 
301 eigenvectors [m+n*number_of .occupations] 

= (m==n) ? cl_f loat ( 1 . , pr ec is ion ) : zero; 

} 

} 

306 /* Execute the "TriDiagonal QL Implicit" Algorithm - from Numerical Receives in C */ 

tql i ( diag , offdiag , number _of _oc cup at ions, eigenvectors, true); 

//Evaluate the partition function 
for(int n=0; n< number _ of _ o c cupat i ons ; n++) { 
311 *partition = *partition + exp(cl_float(diag[n], precision)); 

} 

♦expectation = zero; 

316 //Evaluate the numerator , Tr [ ( a + a*) exp {f (n) + r*b eta * ( a + a*) } ] 

for(int n=0; n< number _ of _ oc cupat i ons ; n++) { 

for(int k = 0; k<number_of .occupations -1 ; k + + ) { //less one 
temp = temp + cl.f loat (2* sqrt (k+1) , precision) 

*cl_float (eigenvectors [k+n* number _ of _ occupations] , precision) 
321 *cl_float (eigenvectors [k+l+n*nMber_of .occupations] , precision); 

} 

♦ expectation = *expectation + t emp * exp ( cl.f loat ( diag [n] , precision)); 
t emp = zer o ; 

} 

326 *expectation = *expectation / *partition ; 

} 



331 /* 

* Find value of partition function for the approximating hamiltonian - hence r dependence. 

* -> could also be known as " cy cl e_ dens i ty _ denominat or " to match below function 
* 

* Evals: Tr [ -\beta( \lambda n(n-l) + (l-\mu)n + r (a+a*) ] 

336 */ 

cl.F part it ion.of _appr ox ( double mu , double r) 
{ 
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cl_F partition = zero; 

//Construct the matrix - first with q=0 
for(int n=0; n< number _ of _ o c cupat i ons ; n++){ 

diag [n] = f_approx(n, 0, mu) ; 

if (n ! =0) //offdiag[0] not used 

offdiag [n] = cl_float (r*beta,precision)*sqrt(cl_float (n, precision)) ; 

} 

/* Execute the "TriDiagonal QL Implicit" Algorithm - from Numerical Receipes in C */ 
tqli (diag , offdiag , number_of .occupations , eigenvectors , false) ; 

//Sum the exponentiated eigenvalues 

for(int n = 0; n< number _ of _ o c cupat i ons ; n + + ) { 

partition = partition + exp(cl_float(diag[n],precision)); 

} 

return partition; 



356 } 



/* 

361 * Evaluate the density of cycles of length q in the approximating hamiltonian case 
* 

* Evals : Tr [ -\beta ( \ lamb da (n + q) (n + q-1) + (l-\mu) (n + q) + r(a + a*)J 
*/ 

cl_F cycle_density_numerator ( int q, double mu , double r) 

366 { 

cl_F numerator = zero ; 

//Construct the matrix - with q non-zero!!! 
for(int n=0; n< number _ of _ oc cupat i ons ; n++) { 
371 diag[n] = f_approx(n, q, mu) ; 

if (n!=0) //offdiag[0] not used 

offdiag [n] = cl_float (r*beta,precision)*sqrt(cl_float (n, precision)) ; 

} 

376 /* Execute the "TriDiagonal QL Implicit" Algorithm - from Numerical Receives in C */ 

tqli (diag, offdiag, number _of _oc cupat i ons , eigenvectors, false); 

//Sum the exponentiated eigenvalues 
for(int n=0; n< number _ of _ o c cupat i ons ; n++) { 
381 numerator = numerator + exp ( cl_f loat ( di ag [n] , pr ec is ion ) ) ; 

} 

return numerator; 



/* 

* Diagonal values of B~H hamiltonian (off-diagonal values generated separately) 

*/ 

cl_F f ( int n, double mu) 

{ 

return cl_f loat (beta , precision )* ( cl_float((mu + lambda - l)*n, precision) 

- cl_f loat ( lambda *n*n , pre ci s ion ) ); 

} 



/* 

* Diagonal values of Approximator of the B-H hamiltonian - with cycles 
401 */ 

cl_F f_approx(int n, int q, double mu) 
{ 

return cl_f loat (-beta , precision )* ( cl_f loat ( lambda *( n+q) *( n+q- 1 ) , precision) 

+ cl_f loat ( ( 1 -mu ) * (n+q) , pr ec i s i on ) ); 

406 } 
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/* 

* Density evaluator 

*/ 

double dens ity ( double r, double mu) 
{ 

cl_F e.vec , temp; 

cl_F numerator = zero ; 

cl_F denominator = zero; 

for(int n=0; n< number _ of _ o c cupat i ons ; n++) { 
diag[n] = f(n, mu) ; 
if (n ! =0) //offdiag[0] not used 

offdiag [n] = cl_float (r*beta,precision)*sqrt(cl_float (n, precision)) ; 
//Create identiy matrix 

for(int m=0; i<nuiber_of .occupations ; m++) { 
eigenvectors [i+n*nuiDber_of .occupations] 

= (m==n) ? cl_f loat ( 1 . , pr ec is ion ) : zero; 

} 

} 

/* Execute the "TriDiagonal QL Implicit" Algorithm - from Numerical Receives in C */ 
//Only want eigenvalues , so set last argument in tqli as true 
tqli(diag, offdiag, number.of .occupations , eigenvectors, true); 

for(int k=0; k< number _ of _ oc cupat i ons ; k++) { 

e_vec = exp ( cl_f loat ( di ag [k] , pre ci sion )) ; 

denominator = denominator + cl.float (e.vec .precision) ; 

t emp = zer o ; 

f or ( int n=l; n<nuiber_of .occupations ; n++) { //n==0 always zero so skip 

temp = temp + n*cl_f loat (eigenvectors [n+k*number_of .occupations] .precision ) 

*cl_f loat (eigenvectors [n+k*number_of .occupations] .precision) ; 

} 

numerator = numerator + cl.float (e_vec*temp , precision ) ; 
return double _appr ox ( numerator /denominator ) ; 



446 

/* Based on the ' tqli ' a I gorithm from Numeric a I Recipes in C 

* Finds the eigenvalues and eigenvectors of a symmetric tridiagonal matrix 

* --enhanced for arbitrary precision arithmetic - reliable up to 10 decimal places it 

appears - - 

*/ 

451 void tql i ( vector < cl.F > fed, vector <cl_F> &e , int n, vector <cl_F> &z , bool e val _ e igenve c s ) 
{ 

int m , 1 , iter , i , k ; 

cl.F s=zero ,r=zero .p=zero .g=zero ,f=zero ,dd=zero ,c=zero ,b=zero; 

456 f or ( i=l ; i <n ; i ++) e[i-l]=e[i]; //shift off-diagonal entries by one 

e [n-1] =zero ; 

f or (1=0 ; Kn ; 1+ + ) { //main loop 
iter =0 ; 
461 do{ 

for(m=l;m<n-l;m++){ 

dd = abs (d [m] ) +abs (d [m + 1] ) ; 

if Cabs (e [m] ) < TINY*dd) break; //this offdiagonal negligible, go onto next 

} 

466 if(m != 1){ 

if (iter++ == 30){ //This algorithm quickly converges , but if not 
verbose << "Too many iterations in TQLI" << endl ; 
throw 1; 

} 

471 if ( e[l] == zero ) throw 1; //Divide-by-zero catch 

g=(d[l+l]-d[l])/(2.0*e[l]) ; 
r=pythag (g , 1 . 0) ; //pythag 
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APPENDIX C. NUMERICAL SOURCE CODE 



g=d[m]-d[l]+e[l]/(g+SIGN(r,g)) ; 
s=c=cl_float (l.precision) ; 
p=zero ; 

for(i=m-l;i>=l;i--){ 
f = s*e [i] ; 
b = c*e [i] ; 

e [i+1] =(r=pythag (f , g) ) ; 
if(r == zer o ) { 

d [i + 1] = d [i + 1] -p ; 

e [m] =zero ; 

break ; 

} 

s = cl_f loat ( f /r , pr ec i s ion ) ; //helps catch underflows 

c=cl_float (g/r , precision ) ; 

g = d [i + 1] -p ; 

r=(d[i]-g)*s+2.0*c*b; 

d [i + 1] =g+(p = s*r) ; 

g=c*r-b ; 

// Next loop can be omitted if eigenvectors not wanted 
if (eval_eigenvecs){ 

for (k=0 ; k<n ; k++) { 

f =z [k+(i+l) *n] ; 

z[k+(i+l)*n]=s*z [k+i*n] +c*f ; 

z[k+i*n]=c*z[k+i*n]-s*f; 

} 

} 

} 

if (r == zero && i >=1) continue; 

d[l]=d[l]-p; 

e[l]=g; 

e [m] =zero ; 



} 

}while(m != 1) 



511 

/* 

* pythag - finds sqrt (a*a + b*b) in best way possible to avoid overflows 

* Needed by tqli() 
*/ 

516 cl_F pythag ( const cl_F a, const cl_F b) 
{ 

cl_F absa=abs(a), absb=abs(b); 

return( (absa > absb) ? absa* sqrt ( 1 . 0+ ( absb * absb ) / ( abs a* absa ) ) : 

(absb==zero ? zero : absb * sqrt ( 1 . 0+ ( absa* absa )/( absb * absb )))) ; 

521 } 
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